---
title: "CE083 - Estatística Computacional I"
author: "Fernando Mayer"
date: "15 de dezembro, 2015"
output: html_document
---

# Sobre R Markdown

- Veja o material disponível em [Rmarkdown-rautu](http://git.leg.ufpr.br/fernandomayer/Rmarkdown-rautu)
- Veja apresentação do [R Users Group](http://www.leg.ufpr.br/~fernandomayer/useR/Rmarkdown-rautu.html)
- Veja o vídeo da apresentação no [Youtube](https://www.youtube.com/watch?v=qFh4XtckP2k)

# Intervalos de confiança

Nesta sessão vamos verificar como utilizar o R para obter intervalos de
confiança para parâmetros de distribuições de probabilidade.

Para fins didáticos, mostrando os recursos do R, vamos mostrar três
possíveis soluções:

- Fazendo as contas passo a passo, utilizando o R como uma calculadora
- Escrevendo uma função
- Usando uma função já existente no R

## Média de uma distribuição normal com variância desconhecida

Considere o seguinte problema:

O tempo de reação de um novo medicamento pode ser considerado como tendo
distribuição Normal e deseja-se fazer inferência sobre a média que é
desconhecida obtendo um intervalo de confiança.  Vinte pacientes foram
sorteados e tiveram seu tempo de reação anotado. Os dados foram os
seguintes (em minutos):

```
2.9  3.4  3.5  4.1  4.6  4.7  4.5  3.8  5.3  4.9
4.8  5.7  5.8  5.0  3.4  5.9  6.3  4.6  5.5  6.2
```

Entramos com os dados com o comando

```{r}
tempo <- c(2.9, 3.4, 3.5, 4.1, 4.6, 4.7, 4.5, 3.8, 5.3, 4.9,
            4.8, 5.7, 5.8, 5.0, 3.4, 5.9, 6.3, 4.6, 5.5, 6.2)
```

Sabemos que o intervalo de confiança para média de uma distribuição
normal com variância desconhecida, para uma amostra de tamanho $n$ é
dado por:

$$
\left(\bar{x} - t_t \sqrt{\frac{s^2}{n}} \quad, \quad
  \bar{x} + t_t \sqrt{\frac{s^2}{n}} \right)
$$

onde $t_t$ é o quantil de ordem $1-\alpha/2$ da distribuição $t$ de 
Student, com $n-1$ graus de liberdade.

Vamos agora obter a resposta das três formas diferentes mencionadas
acima.

### Fazendo as contas passo a passo

Nos comandos a seguir calculamos o tamanho da amostra, a média e a
variância amostral.

```{r}
n <- length(tempo)
n
t.m <- mean(tempo)
t.m
t.v <- var(tempo)
t.v
```

 
A seguir montamos o intervalo utilizando os quantis da distribuição $t$,
para obter um IC a 95\% de confiança. 

```{r}
t.ic <- t.m + qt(c(0.025, 0.975), df = n - 1) * sqrt(t.v/length(tempo))
t.ic
```

### Escrevendo uma função

Podemos generalizar a solução acima agrupando os comandos em uma função.
Nos comandos abaixo, primeiro definimos a função, inclusive adicionando
um argumento `conf` que permitirá que o usuário especifique o nível de
confiança desejado

```{r}
ic.m <- function(x, conf = 0.95){
  n <- length(x)
  media <- mean(x)
  variancia <- var(x)
  tinf <- (1 - conf)/2
  tsup <- 1 - (1-conf)/2
  quantis <- qt(c(tinf, tsup), df = n - 1)
  ic <- media + quantis * sqrt(variancia/n)
  return(ic)
}
```

A seguir, podemos aplicar nossa função ao objeto `tempo`, calculando os
intervalos com 95% e 99% de confiança

```{r}
## IC com 95% de confiança
ic.m(tempo)
## IC com 99% de confiança
ic.m(tempo, conf = 0.99)
```

Escrever uma função é particularmente útil quando um procedimento será
utilizado várias vezes.

### Usando a função `t.test()`

Mostramos as soluções acima para ilustrar a flexibilidade e o uso do
programa. Entretanto não precisamos fazer isto na maioria das vezes
porque o R já vem com várias funções para procedimentos estatísticos
já escritas. Neste caso a função `t.test()` realiza o teste $t$ para a 
média, mas como resultado ela também retorna o intervalo de confiança.
Note que as saídas abaixo coincidem com os obtidos anteriormente.

```{r}
## IC com 95% de confiança
t.test(tempo)
## IC com 99% de confiança
t.test(tempo, conf.level = 0.99)
```

# Testes de hipótese

Nesta sessão vamos verificar como utilizar o R para fazer testes de
hipótese sobre parâmetros de distribuições para as quais os resultados
são bem conhecidos.

Os comandos e cálculos são bastante parecidos com os vistos em
intervalos de confiança e isto nem poderia ser diferente visto que
intervalos de confiança e testes de hipótese são relacionados.

Assim como fizemos com intervalos de confiança, aqui sempre que
possível, e para fins didáticos mostrando os recursos do R, vamos
mostrar três possíveis soluções:

- Fazendo as contas passo a passo, utilizando o R como uma calculadora
- Escrevendo uma função
- Usando uma função já existente no R

## Procedimentos gerais para um teste de hipótese

1. Definir a hipótese nula ($H_0$) e a alternativa ($H_a$)
2. Definir um nível de **significância** $\alpha$ (ex.: $\alpha =
  0,05$), que irá determinar o nível de **confiança** $100(1-\alpha)\%$
  do teste
3. Determinar a **região de rejeição** com base no nível de
  significância $\rightarrow$ valor crítico
4. Calcular a **estatística de teste**, sob a hipótese nula
  $\Rightarrow$ valor calculado 
5. Rejeitar a hipótese nula se a estatística de teste calculada
  estiver dentro da região de rejeição $\Rightarrow$ |valor calculado| >
  |valor crítico|

##  Testes de hipótese para a média: $\sigma$ desconhecido

### Procedimentos gerais para um teste de hipótese com $\sigma$ desconhecido

1. Definir a hipótese nula ($H_0$) e a alternativa ($H_1$)
2. Definir um nível de **significância** $\alpha$ (ex.: $\alpha =
  0,05$), que irá determinar o nível de **confiança** $100(1-\alpha)\%$
  do teste
3. Determinar a **região de rejeição** com base no nível de
  significância $\rightarrow$ valor crítico
4. Calcular a **estatística de teste**, sob a hipótese nula
$$
	t_{calc} = \frac{\bar{x} - \mu_0}{s/\sqrt{n}}
$$
5. Rejeitar a hipótese nula se a estatística de teste calculada
  estiver dentro da região de rejeição ($|t_{calc}| > |t_{crit}|$)

Considerando o problema anterior, do tempo de reação de um medicamento,
faça um teste de hipótse para avaliar se o tempo médio é de 4,5
minutos.

Com isso, já temos as hipótese nula e alternativa para a verdadeira
média $\mu$ do tempo de reação:

$$
\begin{align}
\text{H}_0: \mu = 4,5 \\
\text{H}_1: \mu \neq 4,5 \\
\end{align}
$$

Novamente vamos resolver esse problema de três maneiras diferentes.

### Fazendo as contas passo a passo

Primeiro espacificamos o nível de significância $\alpha = 0,05$, e
achamos o valor crítico de $t$. Depois calculamos a estatística de
teste, e o p-valor associado.

```{r}
## Usando alpha = 0,05, achamos o valor crítico de t com
alpha <- 0.05
tcrit <- abs(qt(alpha/2, df = n - 1))
tcrit
## Estatística de teste
tcalc <- (mean(tempo) - 4.5)/(sd(tempo)/sqrt(n))
tcalc
## tcalc > tcrit?
tcalc > tcrit
## Calculano o p-valor
## Note que essa probabilidade deve ser multiplicada por 2, pois a
## função pt() só calcula uma cauda, e o teste é bicaudal
pt(tcalc, df = n - 1, lower.tail = FALSE) * 2
```

Portanto, por esse resultado, a hipótese nula não pode ser rejeitada.

### Criando uma função

Podemos agregar os comandos acima e algumas definições para criar nossa
própria função para fazer o teste $t$. Siga os passos nos procedimentos
gerais descritos acima e tente criar sua própria função.

```{r, eval=FALSE, echo=FALSE}
## Escreva sua função aqui
teste.t <- function(x, mu, conf = 0.95){
    ...
}
```

```{r}
## Função para fazer o teste t
teste.t <- function(x, mu, conf = 0.95){
    n <- length(x)
    media <- mean(x)
    desvpad <- sd(x)
    tcrit <- abs(qt((1 - conf)/2, df = n - 1))
    tcalc <- (media - mu)/(desvpad/sqrt(n))
    pvalor <- pt(tcalc, df = n - 1, lower.tail = FALSE) * 2
    saida <- list("t critico" = tcrit,
                  "t calculado" = tcalc,
                  "tcalc > tcrit" = tcalc > tcrit,
                  "p-valor" = pvalor)
    return(saida)
}
```

Aplicando a função, temos como resultado

```{r}
teste.t(x = tempo, mu = 4.5)
```

Podemos alterar o nível de significância $\alpha$ para 0,01 por exemplo

```{r}
teste.t(x = tempo, mu = 4.5, conf = 0.99)
```

### Usando a função `t.test()`

Como vimos anteriormente, a função `t.test()` realiza o teste $t$, da
mesma forma que programamos acima. Veja os resultados da saída dessa
função e compare com os que você calculou acima.

```{r}
## Teste t com alpha = 0,05
t.test(tempo, mu = 4.5)
## Teste t com alpha = 0,01
t.test(tempo, mu = 4.5, conf.level = 0.99)
```

## Exercícios

Carregue os dados do endereço abaixo

```{r}
url <- "http://www.leg.ufpr.br/~fernandomayer/dados/crabs.csv"
```

E faça o que se pede abaixo:

1. Testar a hipótese de que a média ($\mu$) de CL é igual a 30 mm (com
  95\% de confiança). As hipóteses são
$$
\begin{align}
H_0: \mu = 30 \\
H_1: \mu \neq 30
\end{align}
$$
2. Testar a hipótese de que a média ($\mu$) de CL é menor ou igual a 30
mm (com 95\% de confiança). As hipóteses são
$$
\begin{align}
H_0: \mu \leq 30 \\
H_1: \mu > 30
\end{align}
$$
3. Testar a hipótese de que a média ($\mu$) de CL é maior ou igual a 30
    mm (com 95\% de confiança). As hipóteses são
$$
\begin{align}
H_0: \mu \geq 30 \\
H_1: \mu < 30
\end{align}
$$
4. Faça dois histogramas lado-a-lado da medida CL para cada uma das
espécies.
5. Testar a hipótese de que a **diferença** entre a média de
CL da espécie azul ($\mu_A$) e a média de CL da espécie laranja
($\mu_L$) é igual a 0 (zero) (com 95\% de confiança). As hipóteses são
$$
\begin{align}
H_0: \mu_A - \mu_L = 0 \quad \Rightarrow \quad \mu_A = \mu_L \\
H_1: \mu_A - \mu_L \neq 0 \quad \Rightarrow \quad \mu_A \neq \mu_L
\end{align}
$$
6. Testar a hipótese de que a **diferença** entre a média de
CL da espécie azul ($\mu_A$) e a média de CL da espécie laranja
($\mu_L$) é **menor** ou iual a 0 (zero) (com 95\% de confiança). Em
outras palavras: *O CL médio é menor para a espécie azul?* As hipóteses
são
$$
\begin{align}
H_0: \mu_A - \mu_L \leq 0 \quad \Rightarrow \quad \mu_A \leq \mu_L \\
H_1: \mu_A - \mu_L > 0 \quad \Rightarrow \quad \mu_A > \mu_L
\end{align}
$$
