Lista 04

library(dplyr)
library(magrittr)
library(readr)
library(ggplot2)

Questão 1

(a)

IC_q1a <- function(x, sigma2, confianca) {
  
  amostra = na.omit(x)
  
  n = length(amostra)
  x_barra = mean(amostra)
  
  alfa = 1-confianca
  z = qnorm(1-alfa)
  
  L = (x_barra - z*sqrt(sigma2/n)) %>% round(1)
  
  cat("IC(m,", confianca*100, "%) = [", L, ",Inf).", sep = "")
  
}

(b)

IC_q1b <- function(x, sigma2, confianca) {
  
  amostra = na.omit(x)
  
  n = length(amostra)
  x_barra = mean(amostra)
  
  alfa = 1-confianca
  z = qnorm(1-alfa)
  
  L = (x_barra + z*sqrt(sigma2/n)) %>% round(1)
  
  cat("IC(m,", confianca*100, "%) = (-Inf,", L, "].", sep = "")
  
}

(c)

IC_q1c <- function(x, sigma, confianca) {
  
  amostra = na.omit(x)
  
  n = length(amostra)
  x_barra = mean(amostra)
  
  alfa = 1-confianca
  z = qnorm(1-alfa/2)
  
  LI = (x_barra - z*sigma/sqrt(n)) %>% round(1)
  LS = (x_barra + z*sigma/sqrt(n)) %>% round(1)
  
  cat("IC(m,", confianca*100, "%) = [", LI, ",", LS, "].", sep = "")
  
}

(d)

IC_q1d <- function(x, sigma, confianca) {
  
  amostra = na.omit(x)
  
  n = length(amostra)
  x_barra = mean(amostra)
  
  alfa = 1-confianca
  z = qnorm(1-alfa/2)
  
  LI = (x_barra - z*sigma/sqrt(n)) %>% round(1)
  LS = (x_barra + z*sigma/sqrt(n)) %>% round(1)
  
  # cat("IC(m,", confianca*100, "%) = [", LI, ",", LS, "].", sep = "")
  
  cont = ifelse(30>=LI & 30<=LS, 1, 0)
  return(cont)
}

amostra = data.frame(
  replicate(1000, rnorm(25, mean = 30, sd = sqrt(25)))
  )

soma = 0
for (i in 1:length(amostra)) {
  ifelse(
    (IC_q1d(amostra[[i]], 5, 0.9)==1), (soma = soma +1), 0
  )
}

# soma
cat(
  (soma/1000)*100, "%", sep = ""
)

## 90.5%

## Como pode ser observado, a porcentagem é aproximadamente 90%, o que indica que em 90% das vezes o intervalo contém o verdadeiro valor de mu

Questão 2

IC_q2 <- function(x, confianca) {
  
  x = na.omit(x)
  
  n = length(x)
  
  p_chapeu = 1/2
  
  alfa = 1 - confianca
  z = qnorm(confianca + alfa/2)
  
  erro = sqrt(0.5*0.5/n)*z
  
  LI = p_chapeu - erro
  LS = p_chapeu + erro
  
  cat("IC(p,", confianca*100, "%) = [", LI, ",", LS, "]", sep= "")
  
}

Questão 3

base_3 <- read_table2("G:/Meu Drive/UFF/Met. Comp. para Est. II/Listas/Lista 04/Base saude.txt", na = "9")

base_3$HIV %<>% factor(labels = c("Não", "Sim"))
base_3$DST %<>% factor(labels = c("Não", "Sim")) 
base_3$Sexo %<>% factor(labels = c("Feminino", "Masculino"))

Como a função da letra d era uma adaptação da c, usei a função criada na letra c.

(a)

n = length(base_3$Peso); n

## [1] 30

## Como n é suficientemente grande, pelo TCL podemos estimar a variância pela variância amostral

IC_q1c(base_3$Peso, sigma = sqrt(var(base_3$Peso, na.rm = T)), confianca = 0.97)

## IC(m,97%) = [138.5,371.1].

## Foi realizada a suposição de que n é suficientemente grande para usar o TCL.

(b)

base_q3b <- base_3 %>% filter(HIV=="Sim") %$% Estatura

## Assumindo que a variância da amostra é um bom estimador para a variância da população

IC_q1c(base_q3b, sigma = sqrt(var(base_q3b, na.rm = T)), confianca = 0.95)

## IC(m,95%) = [76.7,139.5].

## Isso quer dizer que aproximadamente em 95% das vezes o intervalo conterá o verdadeiro valor de mu.

## Foi realizada a suposição de que a variância amostral é um bom estimador para a variância populacional.

(c)

Análogo a b

(d)

base_q3d <- base_3 %>% filter(Sexo=="Feminino", HIV=="Não") %$% Idade

## Assumindo que a variância amostral é um bom estimador para a variância populacional

IC_q1c(base_q3d, sigma = sqrt(var(base_q3d, na.rm = T)), confianca = 0.90)

## IC(m,90%) = [26.4,42.8].

## Isso quer dizer que em 90% das vezes o intervalo conterá o verdadeiro valor de mu.

## Foi realizado a suposição de que a variância amostral é um bom estimador para a variância populacional.

Questão 4

IC_q7 <- function(x, y, conf){
  
  grupos = levels(y)
  
  for(i in 1:length(grupos)){
    
    aux = x[y==grupos[i]] %>% .[!is.na(.)]
    
    cat(grupos[i]," - ", "")
    IC_q1c(aux, sigma = sqrt(var(aux, na.rm = T)), conf)
    cat("\n")
  }
}

IC_q7(base_3$Peso, base_3$Sexo, 0.90)

## Feminino  -  IC(m,90%) = [59.1,312.1].
## Masculino  -  IC(m,90%) = [183.4,423.9].

Questão 5

IC_q7(base_3$Estatura, base_3$DST, conf = 0.95)

## Não  -  IC(m,95%) = [56.2,147.2].
## Sim  -  IC(m,95%) = [97,158.4].

## Para os pacientes observados divididos em grupo possui e não possui DST, pode-se observar o intervalo que conterá o verdadeiro valor da média da estatura de cada grupo com um nível de confiança de 95%.

## Com base nos intervalos dos grupos, leva-nos a supor que a média da Estatura de quem possui DST é maior se comparado com a média de quem não possui DST.

Questão 6

IC_q6_var <- function(x, conf) {
  
  x = na.omit(x)
  
  n = length(x)
  alfa = 1-conf
  
  S2 = var(x, na.rm = T)
  
  x1 = qchisq(alfa/2, df = n-1, lower.tail = F)
  x2 = qchisq(1-alfa/2, df = n-1, lower.tail = F)
  
  LI = (n-1)*S2/x1
  LS = (n-1)*S2/x2
  
  cat(
    "IC(sigma^2,", conf*100, "%) = [", LI %>% round(2), ",", LS %>% round(2), "].", sep = ""
  )
  
}

(a)

base_3 %>% filter(Sexo=="Masculino", HIV=="Sim") %>% select(Idade) %>% pull %>% IC_q6_var(0.90)

## IC(sigma^2,90%) = [28.5,102.72].

## Se realizado o intervalo de confiança para a variância para todas as possíveis amostras da população, teriamos que aproximadamente 90% desses intervalos conteriam o real valor do parâmetro.

## Foram realizadas as suposições de que a amostra veio de uma distribuição normal e que a variância amostral é um bom estimador para a variância populacional.

(b)

IC_q6_var(base_3$Estatura, 0.95)

## IC(sigma^2,95%) = [3076.14,8934.48].

## Se realizado o intervalo de confiança para a variância para todas as possíveis amostras da população, teríamos que aproximadamente 90% desses intervalos conteriam o real valor do parâmetro.

## Foram realizadas as suposições de que a amostra veio de uma distrivui

(c)

base_3 %>% filter(Escol==0)

## # A tibble: 1 x 10
##   Codigo Datacol  Sexo      Idade  Peso Estatura HIV   Escol DST    Tipo
##   <chr>  <chr>    <fct>     <int> <dbl>    <dbl> <fct> <int> <fct> <int>
## 1 AB15   01/01/14 Masculino    18    65      166 Sim       0 Não      NA

## Não é possível encontrar o intervalo de confiança pois o tamanho da amostra é 1.

Questão 7

base_7 <- read_rds("G:/Meu Drive/UFF/Met. Comp. para Est. II/Listas/Lista 04/exames medicos.rds")

(a)

## HDL

## Empiricamente, não é razoável supor normalidade.

ggplot(base_7, aes(sample = base_7$HDL)) + geom_qq() + geom_qq_line()

## LDL

## Empiricamente, é razoável supor normalidade.

ggplot(base_7, aes(sample = base_7$LDL)) + geom_qq() + geom_qq_line()

## Glicose

## Empiricamente, é razoável supor normalidade

ggplot(base_7, aes(sample = base_7$glicose)) + geom_qq() + geom_qq_line()

## Linfocitos

## Empiricamente, não é razoável supor normalidade.

ggplot(base_7, aes(sample = base_7$linfocitos)) + geom_qq() + geom_qq_line()

(b)

IC_media_q7 <- function(x, conf){
  
  x_barra = mean(x, na.rm = T)
  s2 = var(x, na.rm = T)
  
  n = length(x[!is.na(x)])
  alfa = 1-conf
  
  erro = qt(1-alfa/2, df = n-1) * sqrt(s2/n)
  
  LI = x_barra - erro 
  LS = x_barra + erro

  cat(
    "IC (m,", conf*100, "%) = [", LI, ";", LS, "]", sep = ""
  )
}

IC_media_q7(base_7$LDL, 0.98)

## IC (m,98%) = [29.06364;31.40836]

IC_media_q7(base_7$glicose, 0.98)

## IC (m,98%) = [46.36625;51.08975]

(c)

## Empiricamente. Igual feito na letra a.

(d)

IC_var_q7 <- function(x, conf){
  
  n = length(x[!is.na(x)])
  s2 <- var(x, na.rm = T)
  
  alfa = 1-conf
  
  x1 = qchisq(alfa/2, df = n-1, lower.tail = F)
  x2 = qchisq(1-alfa/2, df = n-1, lower.tail = F)
  
  LI = (n-1)*s2 / x1
  LS = (n-1)*s2 / x2
  
  cat(
    "IC (sigma,", conf*100, "% = [", LI, ";", LS, "]", sep = ""
  )
  
}

IC_var_q7(base_7$LDL, 0.95)

## IC (sigma,95% = [3.373065;10.70688]

IC_var_q7(base_7$glicose, 0.95)

## IC (sigma,95% = [13.68889;43.45165]

(e)

IC_prop_q7 <- function(num_s, n, conf, conservador=FALSE){
  if_else(
    conservador==FALSE, p_chapeu = num_s/n, p_chapeu = 1/2
  )
  
  alfa = 1-conf
  
  erro = qnorm(1-alfa/2) * sqrt( (pchapeu*(1-p_chapeu))/n )
  
  LI = p_chapeu - erro
  LS = p_chapeu + erro
  
  cat(
    "IC (p,", conf*100, "%) = [", LI, ";", LS, "]", sep = ""
  )
  
}

QUestão 8

base_8 <- read_rds("G:/Meu Drive/UFF/Met. Comp. para Est. II/Listas/Lista 04/colesterol.rds")

(a)

## Empiricamente, não é razoável supor normalidade.

ggplot(base_8, aes(sample = base_8$HDL)) + geom_qq() + geom_qq_line()

## Empiricamente, é razoável supor normalidade.

ggplot(base_8, aes(sample = base_8$LDL)) + geom_qq() + geom_qq_line()

(b)

## Suponde normalidade independente do resultado na questão acima.

## Supondo normalidade, sim é possível obter o Intervalo de Confiança para a média.

IC_media_q8 <- function(x, conf){
  
  n = length(x[!is.na(x)])
  x_barra = mean(x, na.rm = T)
  s2 = var(x, na.rm = T)
  
  alfa = 1-conf
  
  erro = qt(1-alfa, df = n-1) * sqrt(s2/n)
  
  LI = x_barra - erro
  LS = x_barra + erro
  
  cat(
    "IC (m,", conf*100, ") = [", LI, ";", LS, "]", sep = ""
  )
}

## Isso quer dizer que aproximadamente 95% das vezes o intervalo conterá o verdadeiro valor de mu.

base_8$HDL %>% IC_media_q8(0.95)

## IC (m,95) = [0.2874198;0.3948024]

## Isso quer dizer que aproximadamente 95% das vezes o intervalo conterá o verdadeiro valor de mu.

base_8$LDL %>% IC_media_q8(0.95)

## IC (m,95) = [49.69716;50.44951]

(c)

## Empiricamente, igual feito na letra a.

(d)

## Suponde normalidade independente do resultado na questão acima.

## Supondo normalidade, sim é possível obter o Intervalo de Confiança para a variância.

IC_var_q8 <- function(x, conf){
  
  n = length(x[!is.na(x)])
  s2 = var(x, na.rm = T)
  
  alfa = 1-conf
  
  x1 = qchisq(alfa/2, df = n-1)
  x2 = qchisq(1-alfa/2, df = n-1)

  
  LI = (n-1)*s2/x2
  LS = (n-1)*s2/x1
  
  cat(
    "IC (sigma,", conf*100, ") = [", LI, ";", LS, "]", sep = ""
  )
}

## Isso quer dizer que aproximadamente 95% das vezes o intervalo conterá o verdadeiro valor de sigma.

base_8$HDL %>% IC_var_q8(0.95)

## IC (sigma,95) = [0.07144162;0.128993]

## Isso quer dizer que aproximadamente 95% das vezes o intervalo conterá o verdadeiro valor de sigma.

base_8$LDL %>% IC_var_q8(0.95)

## IC (sigma,95) = [3.50696;6.332068]

Lista 05

library(dplyr)
library(ggplot2)

Questão 1

(a)

## Amostras a. s. obtidas de uma população normal. 
## Variância conhecida.

## H0: mu = mu0
## H1: mu < mu0

TH_media_q1a <- function(x, mu0, sigma2, alfa){
  
  amostra = na.omit(x)
  
  n = length(amostra)
  z = qnorm(1-alfa)
  
  x_barra = mean(amostra)
  x_critico = mu0 - z*sqrt(sigma2/n)
  
  decisao = ifelse((x_barra < x_critico), "Rejeita H0.", "Não rejeita H0.")
  
  cat(
    "Estimativa pontual: ", x_barra, "\n",
    "Região Crítica: [", -Inf, ";", x_critico, "] \n",
    "Decisão: ", decisao, 
    sep = ""
  )
}

## mu0 = 500
rnorm(1000, 500, sqrt(25)) %>% TH_media_q1a(500, 25, 0.05)

## Estimativa pontual: 500.0977
## Região Crítica: [-Inf;499.7399] 
## Decisão: Não rejeita H0.

(b)

## Amostras a. s. obtidas de uma população normal.
## Variância desconhecida.

## H0: sigma2 = sigma0
## H1: sigma2 > sigma0

TH_sigma2_q1b <- function(x, sigma0, alfa){
  
  amostra = na.omit(x)
  n = length(amostra)
  
  s2 = var(amostra)
  x = (n-1)*s2 / sigma0
  
  x_critico = qchisq(alfa, n-1, lower.tail = F)
  
  decisao = ifelse(sigma0 > x_critico, "Rejeita H0.", "Não rejeita H0.")
  
  cat(
    "Estimativa pontual: ", x, "\n",
    "Região Crítica: [", x_critico, ";", +Inf, "] \n",
    "Decisão: ", decisao, 
    sep = ""
  )
}

## sigma0 = 25
rnorm(1000, 500, sqrt(25)) %>% TH_sigma2_q1b(25, 0.05)

## Estimativa pontual: 995.1841
## Região Crítica: [1073.643;Inf] 
## Decisão: Não rejeita H0.

(C)

## Amostras a. s. obtidas de uma população normal.

## H0: 
## H1: 

TH_prop_q1c <- function(x, p0, alfa){
  
  amostra = na.omit(x)
  n = length(amostra)
  
  q = qnorm(1-alfa/2)
  k = sqrt( (p0*(1-p0))/n ) * q
  
  LI = p0 - k
  LS = p0 + k
  
  decisao = ifelse(( p0<LI | p0>LS ), "Rejeita H0.", "Não rejeita H0.")
  
  cat(
    "Estimativa pontual: ", mean(amostra), "\n",
    "Região Crítica: [", LI, ";", LS, "] \n",
    "Decisão: ", decisao, 
    sep = ""
  )
}

c(1,1,1,1,1,1,0,0,0,0,1,0,0,1) %>% TH_prop_q1c(., sum(./length(.)), 0.05)

## Estimativa pontual: 0.5714286
## Região Crítica: [0.3122037;0.8306534] 
## Decisão: Não rejeita H0.

(d)

## H0: mu = mu0
## H1: mu < mu0

TH_media_q1d <- function(x, mu0, sigma2, alfa){
  
  amostra = na.omit(x)
  
  n = length(amostra)
  
  x_barra = mean(amostra)
  
  z = (x_barra-mu0) / sqrt(sigma2/n)
  
  q = qnorm(alfa)
  
  decisao = ifelse((z < q), "Rejeita H0.", "Não rejeita H0.")
  
  cat(
    "Estimativa pontual: ", x_barra, "\n",
    "Estimativa de teste: ", z, "\n",
    "Região Crítica: [", -Inf, ";", q, "] \n",
    "Decisão: ", decisao, 
    sep = ""
  )
}

rnorm(1000, 500, sqrt(25)) %>% TH_media_q1d(500, 25, 0.05)

## Estimativa pontual: 499.8875
## Estimativa de teste: -0.711815
## Região Crítica: [-Inf;-1.644854] 
## Decisão: Não rejeita H0.

(e)

## Erro tipo II: Probabilidade de não rej. H0 dado H0 verdadeiro.

erroII_q1e <- function(n, alfa, sigma2, mu0, mu){
  
  # Região Critica
  x_critico = mu0 - sqrt(sigma2/n)*qnorm(1-alfa)
  
  # Padronizando
  z = (x_critico - mu) / sqrt(sigma2/n)
  
  return( pnorm(z, lower.tail = F) )
}

erroII_q1e(1000, 0.05, 25, 500, 500)

## [1] 0.95

(f)

## O gráfico mostra as probabilidades de se não rejeitar H0 para cada mu diferente.

ggplot(data = data.frame(x = c(0, 60)), aes(x = x)) +
  stat_function(fun = erroII_q1e, args = list(n = 20, alfa = 0.05, sigma2 = 25, mu0 = 30)) +
  geom_vline(xintercept = 30, linetype = "dashed", col = "red") +
  ggtitle("Erro tipo II") + xlab(expression(mu)) + ylab(expression(beta(mu))) 

(g)

## A medida que o tamanho da amostra aumenta os valores a esquerda de mu0 ficam menos prováveis, fazendo a curva crescer quando está mais próximo de mu0

ggplot(data = data.frame(x = c(0, 60)), aes(x = x)) +
  
  stat_function(fun = erroII_q1e, args = list(n = 20 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "red") +
  stat_function(fun = erroII_q1e, args = list(n = 30 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "black") +
  stat_function(fun = erroII_q1e, args = list(n = 40 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "purple") +
  stat_function(fun = erroII_q1e, args = list(n = 50 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "green") +
  stat_function(fun = erroII_q1e, args = list(n = 100, alfa = 0.05, sigma2 = 25, mu0 = 30), color = "blue") +
  
  geom_vline(xintercept = 30, linetype = "dashed", col ="red") +
  ggtitle("Erro tipo II") + xlab(expression(mu)) + ylab(expression(beta(mu))) 

(h)

poder_q1h <- function(n, alfa, sigma2, mu0, mu){
  
  # Região Crítica
  q = qnorm(1-alfa/2)
  LI = mu0 - sqrt(sigma2/n)*q
  LS = mu0 + sqrt(sigma2/n)*q
  
  # Padronizando
  z1 = (LI - mu) / sqrt(sigma2/n)
  z2 = (LS - mu) / sqrt(sigma2/n)
  
  return(pnorm(z1) + pnorm(z2, lower.tail = F))
  
}

rnorm(1000, 500, sqrt(25)) %>% length() %>% poder_q1h(0.05, 25, 500, 500)

## [1] 0.05

(i)

## O gráfico mostra as probabilidades de se rejeitar H0 para cada mu diferente.

ggplot(data = data.frame(x = c(0, 60)), aes(x = x)) +
  stat_function(fun = poder_q1h, args = list(n = 20, alfa = 0.05, sigma2 = 25, mu0 = 30)) +
  geom_vline(xintercept = 30, linetype = "dashed", col ="red") +
  ggtitle("Função Poder") + xlab(expression(mu)) + ylab(expression(pi(mu))) 

(j)

## A medida que o tamanho da amostra aumenta maior é o poder do teste, ou seja, a probabilidade de se rejeitar valores diferentes de mu0 aumenta.

ggplot(data = data.frame(x = c(0, 60)), aes(x = x)) +
  
  stat_function(fun = poder_q1h, args = list(n = 20 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "red") +
  stat_function(fun = poder_q1h, args = list(n = 30 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "black") +
  stat_function(fun = poder_q1h, args = list(n = 40 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "purple") +
  stat_function(fun = poder_q1h, args = list(n = 50 , alfa = 0.05, sigma2 = 25, mu0 = 30), color = "green") +
  stat_function(fun = poder_q1h, args = list(n = 100, alfa = 0.05, sigma2 = 25, mu0 = 30), color = "blue") +
  
  geom_vline(xintercept = 30, linetype = "dashed", col ="red") +
  ggtitle("Função Poder") + xlab(expression(mu)) + ylab(expression(pi(mu))) 

Questão 2

(a)

## Amostras a. s. obtidas de uma população normal com média mu e variância sigma2 desconhecida.

## Suponha que n é grande.

## H0: mu = mu0
## H1: mu > mu0

TH_media_q2a <- function(x, alfa, mu0){
  
  amostra = na.omit(x)
  n = length(amostra)
  
  ## x_barra segue distribuição normal pois é combinação linear de distribuições normais.
  x_barra = mean(amostra)
  
  s2 = var(amostra)
  
  q = qnorm(alfa, lower.tail = F)
  
  k = sqrt(s2/n) * q
  
  x_critico = mu0 + k
  
  decisao = if_else((x_barra > x_critico), "Rejeita H0", "Não rejeita H0")
  
  cat(
    "Estimativa pontual: ", x_barra, "\n",
    "Região Crítica: [", x_critico, " ; ", +Inf, "] \n",
    "Decisão: ", decisao, " a um nível de significância ", alfa,
    sep = ""
  )
}

rnorm(1000, 500, sqrt(25)) %>% TH_media_q2a(0.05, 500)

## Estimativa pontual: 499.989
## Região Crítica: [500.2685 ; Inf] 
## Decisão: Não rejeita H0 a um nível de significância 0.05

(b)

## Amostras a. s. obtidas de uma população com distribuição normal com média mu e variância desconhecida.

## H0: mu = mu0
## H1: mu > mu0

TH_media_q2b <- function(x, alfa, mu0){
  
  amostra = na.omit(x)
  n = length(amostra)
  
  ## x_barra segue distribuição normal pois é combinação linear de distribuições normais.
  x_barra = mean(amostra)

  s2 = var(amostra)
  
  t = (x_barra - mu0) / sqrt(s2/n)
  
  x_critico = qnorm(alfa, lower.tail = F)
  
  decisao = if_else((t > x_critico), "Rejeita H0", "Não rejeita H0")
  
  cat(
    "Estatística de teste: ", t, "\n",
    "Região Crítica: [", x_critico, " ; ", +Inf, "] \n",
    "Decisão: ", decisao, " a um nível de significância ", alfa,
    sep = ""
  )
}

rnorm(1000, 500, sqrt(25)) %>% TH_media_q2a(0.05, 500)

## Estimativa pontual: 499.9543
## Região Crítica: [500.2576 ; Inf] 
## Decisão: Não rejeita H0 a um nível de significância 0.05

(c)

erro_tipoII_q2c <- function(n, alfa, s2, mu0, mu){
  
  ## Região critica 
  
  k = sqrt(s2/n)*qnorm(1-alfa)
  x_critico = mu0 + k
  
  ## Padronizando e calculando probabilidade
  
  z = (x_critico - mu) / sqrt(s2/n)
  
  return( pnorm(z) )
}

erro_tipoII_q2c(1000, 0.05, 25, 500, 500)

## [1] 0.95

(d)

## O gráfico mostra as probabilidades de se não rejeitar H0 para cada mu diferente.

ggplot(data = data.frame(x = c(0, 20)), aes(x = x)) +
  stat_function(fun = erro_tipoII_q2c, args = list(n = 20, alfa = 0.05, s2 = 25, mu0 = 10)) +
  geom_vline(xintercept = 10, linetype = "dashed", col = "red") +
  ggtitle("Erro tipo II") + xlab(expression(mu)) + ylab(expression(beta(mu)))

(e)

## Pode-se afirmar que a medida que o tamanho da amostra aumenta os valores a direita de mu0 ficam menos prováveis.

ggplot(data = data.frame(x = c(0, 60)), aes(x = x)) +
  
  stat_function(fun = erro_tipoII_q2c, args = list(n = 20 , alfa = 0.05, s2 = 25, mu0 = 30), color = "red") +
  stat_function(fun = erro_tipoII_q2c, args = list(n = 30 , alfa = 0.05, s2 = 25, mu0 = 30), color = "black") +
  stat_function(fun = erro_tipoII_q2c, args = list(n = 50 , alfa = 0.05, s2 = 25, mu0 = 30), color = "green") +
  stat_function(fun = erro_tipoII_q2c, args = list(n = 100, alfa = 0.05, s2 = 25, mu0 = 30), color = "blue") +
  
  geom_vline(xintercept = 30, linetype = "dashed", col = "red") +
  ggtitle("Erro tipo II") + xlab(expression(mu)) + ylab(expression(beta(mu)))

(f)

ggplot(data = data.frame(x = c(0, 60)), aes(x = x)) +
  
  stat_function(fun = erro_tipoII_q2c, args = list(n = 20, alfa = 0.01, s2 = 25, mu0 = 30), color = "red") +
  stat_function(fun = erro_tipoII_q2c, args = list(n = 20, alfa = 0.05, s2 = 25, mu0 = 30), color = "black") +
  stat_function(fun = erro_tipoII_q2c, args = list(n = 20, alfa = 0.07, s2 = 25, mu0 = 30), color = "green") +
  stat_function(fun = erro_tipoII_q2c, args = list(n = 20, alfa = 0.10, s2 = 25, mu0 = 30), color = "blue") +
  
  geom_vline(xintercept = 30, linetype = "dashed", col = "red") +
  ggtitle("Erro tipo II") + xlab(expression(mu)) + ylab(expression(beta(mu)))

(g)

poder_q2 <- function(n, alfa, s2, mu0, mu){
  
  ## Região Crítica
  
  k = sqrt(s2/n)*qnorm(1-alfa)
  x_critico = mu0 + k
  
  ## Padronizando e calculando probabilidade
  
  z = (x_critico - mu) / sqrt(s2/n)
  
  return(pnorm(z, lower.tail = F))
}

poder_q2(1000, 0.05, 25, 500, 500)

## [1] 0.05

## Plotando a função 
ggplot(data = data.frame(x = c(0,20)), aes(x = x)) +
  stat_function(fun = poder_q2, args = list(n = 20, alfa = 0.05, s2 = 25, mu0 = 10)) +
  geom_vline(xintercept = 10, linetype = "dashed", col = "red") +
  ggtitle("Erro tipo II") + xlab(expression(mu)) + ylab(expression(pi(mu)))

Questão 3

base_3 <- readRDS("G:/Meu Drive/UFF/Met. Comp. para Est. II/Listas/Lista 05/BaseGenero.rds")

(a)

## H0: mu = 70
## H1: mu > 70

## Sim, com base no nível de significância de 5%, há evidências de que o pesquisador está certo. 

base_3 %>% filter(Sexo=="Homem") %>% select(Peso) %>% 
  t.test(alternativa = "greater", mu = 70) 

## 
##  One Sample t-test
## 
## data:  .
## t = 13.065, df = 24, p-value = 2.107e-12
## alternative hypothesis: true mean is not equal to 70
## 95 percent confidence interval:
##  79.31101 82.80456
## sample estimates:
## mean of x 
##  81.05778

(b)

## H0: mu = 30
## H1: mu < 30

## Não, com base no nível de significância de 1% há evidências de que o pesquisador está errado.

TH <- function(x, mu0, alfa){
  
  amostra = na.omit(x)
  n = length(amostra)
  
  ## Supondo que n é grande o suficiente
  s2 = var(amostra)
  
  k = sqrt(s2/n) * qt(alfa, n-1, lower.tail = F)
  
  x_critico = mu0 - k
  
  cat(
    "RC = [", -Inf, " ; ", x_critico, "]", sep = ""
  )
}

base_3 %>% filter(Sexo=="Homem") %>% select(Idade) %>% pull %>% TH(30, 0.01)

## RC = [-Inf ; 27.73446]

(c)

## H0: mu  =  20
## H1: mu =/= 20

## Rejeita h0; Sim, com base no nível de significância de 3% há evidências de que o pesquisador está correto.

base_3 %>% filter(Sexo=="Mulher") %>% select(Idade) %>% t.test(conf.level = 0.99, alternativa = "two.sided", mu = 20)

## 
##  One Sample t-test
## 
## data:  .
## t = 6.5822, df = 19, p-value = 2.667e-06
## alternative hypothesis: true mean is not equal to 20
## 99 percent confidence interval:
##  22.96809 27.53191
## sample estimates:
## mean of x 
##     25.25

(d)

## H0: p = 65
## H1: p > 65

## Com base nos dados e adotando um nível de significância de 5%, não rejeitamos H0, ou seja, há evidência de que a proporção de mulheres com peso superior a 70kg não é superior a 65%.

amos <- base_3 %>% filter(Sexo=="Mulher") %>% select(Peso) %>% pull
proporcao = ifelse(amos>70, 1, 0)

prop.test(sum(proporcao), length(proporcao), alternative = "greater", p = 0.65, correct = F)

## 
##  1-sample proportions test without continuity correction
## 
## data:  sum(proporcao) out of length(proporcao), null probability 0.65
## X-squared = 14.066, df = 1, p-value = 0.9999
## alternative hypothesis: true p is greater than 0.65
## 95 percent confidence interval:
##  0.1273772 1.0000000
## sample estimates:
##    p 
## 0.25