Análise exploratória e modelação com r parte 3

Minicurso I

Análise Exploratória e Modelação com R 3. Tópicos de Inferência Estatística

Amílcar Oliveira ; Teresa Oliveira

[email protected] [email protected]

DCeT, Universidade Aberta & CEAUL, Universidade de Lisboa

Manaus, 12 a 15 de agosto, 2014

mailto:[email protected]

mailto:[email protected]

Tópicos

1. Distribuições de Probabilidade

2. Intervalos de Confiança

3. Testes de hipóteses

4. Métodos não Paramétricos

2 MINICURSO CIAEEAR

MANAUS (AM), 12-15 AGOSTO 2014

Distribuições de Probabilidade

3

DISTRIBUIÇÕES CONTÍNUAS

DISTRIBUIÇÕES DISCRETAS

MINICURSO CIAEEAR


Distribuição Descrição

rnorm(n,mean,sd)

runif(x,min,max)

rt(x,df)

rf(x)

rexp(x)

rchisq(x,df)

rbinom(x,n,p)

rgeom(x,p)

rpois(x,lambda)

rhyper(x,m,n,k)

Normal Uniforme t de Student F Exponencial Qui-Quadrado Binomial Geométrica Poisson Hipergeométrica


4

DISTRIBUIÇÃO NORMAL

Como gerar valores de uma distribuição normal de média µ e variância 𝜎2 ?

Exemplo: gerar 20 valores de uma distribuição normal de média µ=10 e variância=4

X <- rnorm(20,mean=10,sd=2)

[1] 10.454307 10.109470 10.819990 7.593200 14.135418 6.808652

9.907171 12.377358 9.653927 10.555182 9.588058 12.041869 8.157866

8.119346 8.003958 10.399826 11.971313 13.579755 12.913764

7.725487

O R aceita usar uma função dentro de outra função, assim podemos também resolver da

forma

x <- rnorm(20,mean=10,sd=sqrt(4))

Como todas as funções em geral, a função rnorm() considera valores dos parâmetros

“por defeito”:

x<-rnorm(20) # apresenta 20 valores de uma distribuição normal de

média 0 e variância 1

MINICURSO CIAEEAR



5


Vamos verificar a ideia anterior, gerando 10000 valores de uma Normal (0,1) e calculando

de seguida a média e a variância dos valores gerados.

> x<-rnorm(10000)

> mean(x)

[1] 0.0007968414

> sd(x)

[1] 1.003954

O resultado não dá os valores exatos, pois trata-se da geração de números pseudo-

aleatórios. No entanto a medida que se aumente o número de valores gerados a precisão

aumenta.

MINICURSO CIAEEAR



6


Podemos visualizar a distribuição gerada inicialmente usando a função plot().

> x<-rnorm(20,10,2)

> plot(x)

MINICURSO CIAEEAR


5 10 15 20

68

10

12

Index

x


7


E na forma de histograma usando a função hist().

> hist(x)

MINICURSO CIAEEAR


Histogram of x

x

Fre

qu

en

cy

6 8 10 12 14

01

23

45


8


Usando o gráfico de estremos e quartis (boxplot) boxplot().

> boxplot(x)

máximo

3º quartil

mediana

1º quartil

mínimo

MINICURSO CIAEEAR


68

10

12


9

DISCRETAS

Comecemos por considerar uma função de probabilidade

𝑃 𝑋 = 𝑘 =

0.2 𝑠𝑒 𝑥 = 10.1 𝑠𝑒 𝑥 = 20.4 𝑠𝑒 𝑥 = 30.3 𝑠𝑒 𝑥 = 4

0 𝑜𝑢𝑡𝑟𝑜𝑠 𝑐𝑎𝑠𝑜𝑠

Como representar graficamente no R?

> x<-0:4

> y<-c(0,0.2,0.1,0.4,0.3)

> plot(x,y, type="h")

> points(1,0.2)

> points(2,0.1)

> points(3,0.4)

> points(4,0.3)

MINICURSO CIAEEAR


0 1 2 3 4

0.0

0.1

0.2

0.3

0.4

x

y


10

DISCRETAS

E agora a função acumulada

𝑃 𝑋 ≤ 𝑘 = 𝐹 𝑋 =

0.2 𝑠𝑒 𝑥 < 10.3 𝑠𝑒 𝑥 < 20.7 𝑠𝑒 𝑥 < 31 𝑠𝑒 𝑥 ≥ 4

Como representar graficamente no R?

x<-0:4

y<-c(0,0.2,0.3,0.7,1)

plot(x,y, type="s")

MINICURSO CIAEEAR


0 1 2 3 4

0.0

0.2

0.4

0.6

0.8

1.0

x

y


11

DISTRIBUIÇÃO BINOMIAL

A probabilidade de serem observados k sucessos é dada por

𝑃 𝑋 = 𝑘 =𝑛𝑘

𝑝𝑘(1 − 𝑝)𝑛−𝑘

No R a função que permite gerar números com distribuição binomial é rbinom()

Exemplo: Gerar 50 valores de uma distribuição binomial cuja probabilidade de sucesso é

de 10%, sendo o número de experiências (observações), 150.

> x <- rbinom(50, size=150, p=0.1)

> x

[1] 20 11 16 18 9 12 7 15 8 14 25 15 16 22 12 18 18 18 8 7

15 8 20 17 19

[26] 14 16 11 17 14 15 15 9 13 18 18 17 18 17 13 15 15 13 20 15

20 12 18 18 15

Este resultado poderia ser o de um estudo de epidemiologia, com

150 participantes, com uma prevalência de uma determinada doença

de 10% (correspondendo a uma probabilidade de sucesso de 0.1).

MINICURSO CIAEEAR



12

DISTRIBUIÇÃO BINOMIAL (abordagem gráfica)

Representando através do diagrama de caule e folhas com a função stem()

> stem(x)

The decimal point is at the |

6 | 00

8 | 00000

10 | 00

12 | 000000

14 | 000000000000

16 | 0000000

18 | 0000000000

20 | 0000

22 | 0

24 | 0

MINICURSO CIAEEAR



13


E se quisermos duplicar a escala do gráfico… por exemplo 2 vezes…

> stem(x,scale=2)

The decimal point is at the |

7 | 00

8 | 000

9 | 00

10 |

11 | 00

12 | 000

13 | 000

14 | 000

15 | 000000000

16 | 000

17 | 0000

18 | 000000000

19 | 0

20 | 0000

21 |

22 | 0

23 |

24 |

25 | 0 MINICURSO CIAEEAR



14


Agrupando gráficos …

> par(mfrow=c(2,3))

> hist(rbinom(50, size=200, p=0.02))






Desta forma podemos ver no conjunto, como se altera a distribuição do número de pessoas atingidas pela doença ao aumentarmos o valor da prevalência da doença.

MINICURSO CIAEEAR


Histogram of rbinom(50, size = 200, p = 0.02)

rbinom(50, size = 200, p = 0.02)

Fre

quen

cy

1 2 3 4 5 6 7 8

02

46

810

12


rbinom(50, size = 200, p = 0.04)

Fre

quen

cy

2 4 6 8 10 14

05

1015


rbinom(50, size = 200, p = 0.06)

Fre

quen

cy

6 8 10 14 18

05

1015

20


rbinom(50, size = 200, p = 0.08)

Fre

quen

cy

10 15 20 25

02

46

810


rbinom(50, size = 200, p = 0.1)

Fre

quen

cy

10 15 20 25 30

02

46

810


rbinom(50, size = 200, p = 0.5)

Fre

quen

cy

80 90 100 110

05

1015

Inferência Estatística na comparação de amostras (alguns exemplos)

15

O R possui funções e permite procedimentos para a realização de um vasto número de

testes estatísticos.

Vamos a partir de agora abordar alguns deles e nos quais usaremos os seguintes

conjuntos de dados:

x<-c(16.1,15.2,18.7,10.2,13.6,14.2,19.3,17.7,12.1,13.2,14.2,15.2,15.8,16.0,15.0)

y<-c(13.1,13.6,14.7,12.2,14.2,13.9,15.3,13.7,13.5,14.0,16.9,16.0,13.1,17.0,11.9)

MINICURSO CIAEEAR


Intervalos de Confiança / Testes de hipóteses

16

Teste t (de Student) Para uma média simples

Vamos testar se a média da população (de onde a amostra x foi retirada) é igual a 14. Para realizarmos este teste no R basta usar o comando:

> t.test(x, mu=14, conf.level=0.99)

One Sample t-test

data: x

t = 1.7802, df = 14, p-value = 0.09675

alternative hypothesis: true mean is not equal to 14

99 percent confidence interval:

13.26059 16.93941

sample estimates:

mean of x

15.1

MINICURSO CIAEEAR



17

Interpretação do output

One Sample t-test

data: x

t = 1.7802, df = 14, p-value = 0.09675

alternative hypothesis: true mean is not equal to 14


13.26059 16.93941

sample estimates:

mean of x

15.1

Agora basta fazer a interpretação correta do output do R.

Para saber se a hipótese nula não é rejeitada, basta verificar o valor do p-value

e compará-lo com o nível de significância estipulado, neste caso 1%.

Dado que o valor do p-value é maior que o nível de significância estabelecido,

então não rejeitamos a hipótese nula, o que significa que a população de onde

foi retirada a amostra x possui média estatísticamente igual ao valor 14.

Podemos concluir também que o I.C. para a média populacional, considerando

o mesmo nível de significância de 1% é [13.26059 , 16.93941].

MINICURSO CIAEEAR


18

Teste t (de Student) – duas amostras independentes

Para a realização do teste t, teremos neste caso de admitir ainda dois

pressupostos: populações homocedásticas, e que seguem uma

distribuição normal. Vejamos o seguinte exemplo:

Suponha que queremos testar se “x” e “y” possuem médias

estatisticamente iguais a um nível de significância de 5%. Suponha

também que essas amostras são independentes.

> t.test(x,y,conf.level=0.95)

Welch Two Sample t-test

data: x and y

t = 1.2181, df = 23.806, p-value = 0.2351

alternative hypothesis: true difference in means is not equal to 0


-0.6209537 2.4076204

sample estimates:

mean of x mean of y

15.10000 14.20667 MINICURSO CIAEEAR



19

Interpretação do Output

Dado que o valor do p-value é maior que o nível de significância estabelecido

(5%), então não rejeitamos a hipótese nula, o que significa que estatísticamente

ambas as populações têm a mesma média.

Podemos concluir também que o I.C. para a diferença entre médias populacionais, considerando o mesmo nível de significância de 5% é

[-0.6209537 , 2.4076204]

MINICURSO CIAEEAR


20

Teste t (de Student) – duas amostras dependentes Para a realização do teste t, teremos neste caso de admitir dois

pressupostos: amostras homocedásticas, seguem uma distribuição normal.

Considerando o mesmo exemplo:

Suponhamos que queremos testar se “x” e “y” possuem médias

estatisticamente iguais a um nível de significância de 5%. > t.test(x,y,conf.level=0.95, paired=T)

Paired t-test

data: x and y

t = 1.4431, df = 14, p-value = 0.171

alternative hypothesis: true difference in means is not equal

to 0


-0.4343352 2.2210019

sample estimates:

mean of the differences

0.8933333 MINICURSO CIAEEAR



21

Interpretação do Output

Dado que o valor do p-value é maior que o nível de significância estabelecido

(5%), então não rejeitamos a hipótese nula, o que significa que estatísticamente

ambas as populações têm a mesma média.

Podemos concluir também que o I.C. para a diferença entre médias populacionais, considerando o mesmo nível de significância de 5% é

[-0.4343352 , 2.2210019]

Através deste resultado, podemos também chegar à mesma conclusão, ou seja

estando o valor da diferença entre médias a testar(neste caso é zero) contido

no intervalo leva-nos a concluir pela não rejeição da hipótese nula.

MINICURSO CIAEEAR


22

Teste F > x<-

c(16.1,15.2,18.7,10.2,13.6,14.2,19.3,17.7,12.1,13.2,14.2,15.2,15.

8,16.0,15.0)

> y<-

c(13.1,13.6,14.7,12.2,14.2,13.9,15.3,13.7,13.5,14.0,16.9,16.0,13.

1,17.0,11.9)

> var.test(x,y)

F test to compare two variances

data: x and y

F = 2.4468, num df = 14, denom df = 14, p-value = 0.1055

alternative hypothesis: true ratio of variances is not equal to 1


0.8214631 7.2880075

sample estimates:

ratio of variances

2.4468

MINICURSO CIAEEAR



Métodos não paramétricos

23

As técnicas estatísticas não paramétricas procuram providenciar

procedimentos estatísticos inferenciais menos exigentes do que no caso

paramétrico, ou seja com pressupostos relaxados relativamente às

distribuições subjacentes às populações.

São técnicas de um modo geral bastante simples e portanto acessíveis para

um elevado número de pessoas, não necessariamente estatísticos,

encontrando-se muitas aplicações em áreas como as Ciências Sociais, os

Recursos Humanos e a Psicologia.

Sabe-se que, quanto maior for a informação disponível acerca da distribuição

subjacente aos dados, melhores serão os procedimentos inferenciais

respectivos. Verifica-se que um procedimento não paramétrico geralmente

conduz a uma maior variância (no caso da estimação pontual) e a uma menor

potência (nos testes de hipóteses) quando comparado com um paramétrico.

No entanto muitas são as situações em que é desejável dispensar

pressupostos, sendo nesta altura bastante importante o recurso aos métodos

não paramétricos. MINICURSO CIAEEAR


Métodos não paramétricos

24

De um modo geral os testes não paramétricos são de simples recursos

computacionais, apresentando ainda a vantagem de se poderem aplicar a

situações de variáveis com tipologias pobres em termos de quantificação

ou de ordenação, sendo muitos deles aplicáveis a dados de caracter

qualitativo, ou que apenas se apresentam de acordo com uma escala

ordinal.

Serão abordados os seguintes tipos de testes:

• Teste para a independência (Qui-quadrado)

• Testes para a comparação de médias (O teste U de Mann-Whitney:

comparação de duas amostras independentes ; Kruskal-Wallis

comparação de k amostras independentes)

• Testes de ajustamento (Qui-quadrado, Kolmogorov-Smirnov, Shapiro-

Wilk

• Teste de aleatoriedade (Qui-quadrado)

MINICURSO CIAEEAR


25

Exemplo(teste para a independência):

Num gabinete de planeamento familiar, foi efectuado um estudo, para averiguar

se existe associação entre o estado civil do casal (casados, solteiros) e a

existência de filhos no agregado familiar (com filhos, sem filhos).

Para tal, durante um certo período de tempo, foram registadas observações, que

se apresentam na tabela

Com base nestes valores verifique se de facto existe associação.

Com filhos Sem filhos

Casados 35 25

Solteiros 20 45

MINICURSO CIAEEAR


26

R: Seja X a variável aleatória representativa do estado civil e Y a variável

aleatória representativa da situação (ter ou não ter filhos). As hipóteses a

testar serão H0 : Há independência entre X e Y

H1 : Não há independência entre X e Y

Apresenta-se de seguida o procedimento a seguir no R para resolução deste

exercício:

> count<-matrix(c(35,20,25,45),nrow=2)

> count

[,1] [,2]

[1,] 35 25

[2,] 20 45

> chisq.test(count,correct=F)

Pearson's Chi-squared test

data: count

X-squared = 9.6206, df = 1, p-value = 0.001924

Conclusão: Uma vez que o valor obtido p-value = 0.001924 < 0.05

então rejeita-se H0 , ao nível de significância de 5%. Concluindo-se que não

há independência entre o estado civil e o ter ou não ter filhos

MINICURSO CIAEEAR


27

Exemplo 2 (teste para a comparação de médias):

Na tabela encontram-se os valores das resistências à compressão (MPa) de

provetes cúbicos de betão, obtidos por dois processos distintos de

produção. Com base nesses dois processos retiram-se de forma aleatória,

duas amostras, uma com 12 provetes e outra com 8 provetes. Pretende-se

saber se existe ou não, diferenças significativas entre os dois processos de

produção. Para tal use o teste U de Mann-Whitney, considerando um nível

de significância de 1%.

Amostra A 24.3 23.5 26.4 25.5 27.4 25.8 23.4 25.9 28.5 25.5 29.4 27.2

Amostra B 24.1 25.5 24.8 23.3 24.5 26.7 25.7 27.9

R: Uma das possibilidades para resolver é utilizar o teste de Wilcoxon Rank

Sum Test implementado no R e que é equivalente ao teste U de Mann-

Whitney.

MINICURSO CIAEEAR


28

Assim teremos

> amostraA <-

c(24.3,23.5,26.4,25.5,27.4,25.8,23.4,25.9,28.5,25.5,29.4,27.2)

> amostraB <- c(24.1,25.5,24.8,23.3,24.5,26.7,25.7,27.9)

wilcox.test(amostraA,amostraB)

Wilcoxon rank sum test with continuity correction

data: amostraA and amostraB

W = 61, p-value = 0.3341

alternative hypothesis: true location shift is not equal to 0

Como p-value=0.3341>0.01, não rejeitamos a hipótese nula, ou seja ao

nível de significância considerado (1%), concluímos que não existem

diferenças entre os diferentes sistemas de produção.

MINICURSO CIAEEAR


29

Exemplo 3 (teste de ajustamento):

Escolheram-se 80 pessoas de forma aleatória e registou-se o peso de cada

uma (em quilogramas).

Pretendemos saber se esse conjunto de dados segue a distribuição de qui-

quadrado, considerando-se um nível de significância de 5%.

> pesos

[1] 46.88 47.17 64.46 67.84 85.76 65.41 60.10 75.84 61.21

61.65 63.87 53.95

[13] 63.66 69.06 76.41 75.56 69.04 35.18 66.42 58.78 73.02

51.69 90.88 53.01

[25] 64.31 61.91 79.42 57.78 62.73 60.63 63.29 46.53 84.64

61.76 85.08 59.66

[37] 54.89 94.18 59.89 68.56 75.66 72.06 62.00 43.43 73.38

73.31 66.37 73.72

[49] 66.15 67.79 47.17 64.46 67.84 51.23 34.33 65.43 76.34

87.23 98.23 76.34

[61] 86.34 65.44 45.44 67.99 67.44 65.87 86.29 76.88 55.44

34.67 46.78 34.99

[73] 76.44 56.78 65.87 57.29 73.27 32.44 67.45 68.45

MINICURSO CIAEEAR


30

Agora basta entrar com o comando corretamente:

ks.test(pesos, #amostra a ser testada

“pchisq”, #”p” seguido do nome da distribuição

79) #são os graus de liberdade da amostra

One-sample Kolmogorov-Smirnov test

data: pesos

D = 0.6402, p-value < 2.2e-16

alternative hypothesis: two.sided

Como o p-value é menor ou igual a 0,05 (5%) podemos assumir que os

dados não tem distribuição qui-quadrado ao nível de significância de 5%.

MINICURSO CIAEEAR


31

Exemplo 4 (teste para verificação da normalidade):

> shapiro.test(pesos)

Shapiro-Wilk normality test

data: pesos

W = 0.9777, p-value = 0.1756

Aplicando o teste de Shapiro-Wilk ao exemplo anterior obtemos:

Concluimos (ao nível de significância de 5%) pela não rejeição da hipótese

de normalidade das observações.

MINICURSO CIAEEAR


32

> qqnorm(pesos) #obtendo o normal probability plot só para

comparação

>qqline(pesos) # colocando uma linha auxiliar

-2 -1 0 1 2

30

40

50

60

70

80

90

10

0

Normal Q-Q Plot

Theoretical Quantiles

Sa

mp

le Q

ua

ntile

s

MINICURSO CIAEEAR


33

Exemplo 5 (teste de aleatoriedade):

Suponha que você pretende averiguar se uma amostra de tamanho (n= 200)

com números de 0 a 10 é realmente aleatória.

Para gerar a amostra basta fazer:

amostra<-round(runif(200)*10)

Para obter as freqüências basta usar a função table: >freq<-table(amostra)

>freq

>amostra

0 1 2 3 4 5 6 7 8 9 10

10 18 24 32 17 18 20 18 18 16 9

Uma forma prática de se obter essas frequências é usar o comando

as.numeric().

>chisq.test(as.numeric(freq))

Chi-squared test for given probabilities

data: as.numeric(freq)

X-squared = 21.21, df = 10, p-value = 0.01968 MINICURSO CIAEEAR


Software

Análise exploratória e modelação com r parte 3