#===============================================================================
#. Título: Inferência Estatística - Entender a Origem distribuição t de Student
#. Autor : José Cláudio Faria/UESC/DCET
# Data  : 30/08/2021 18:59:59
#===============================================================================
# Objetivos
#===============================================================================
# a) Gerar uma população normal (parâmetros definidos pelo usuário);
# b) Amostrar repetidamente a população (definido pelo usuário);
# c) Estimar o valor t a partir de cada amostra;
# d) Plotar o histograma dos valores tcal, que deve ser semelhante a função t
#    (parâmetros definidos pelo usuário) se o número de
#    repetições do processo for suficiente para se aproximar do limite;
# e) Sobrepor ao histograma a função densidade de probabilidade específica.
#===============================================================================

#.. ___Ini opções___
mu   <-    40  # M POP
s    <-     5  # DP POP

n    <-    5  # Tamanho da amostra
nsim <-   1e4  # Número de simulações

p    <-  5/100 # Probabilidade adotada na inferência
k    <-    100 # Número de intervalos do histograma
#.. ___Fim opções___

#.. Gerar população e bootstrap para calcular o vetor dos valores t
m <- matrix(rnorm(n * nsim,
                  mu,
                  s),
            nrow=nsim,
            ncol=n)

tm <- apply(m,
            1,
            function(x) (mean(x) - mu)/(sd(x)/sqrt(n)))

#tm <- apply(m,
#            1,
#            function(x) (mean(x) - mu)/(s/sqrt(n)))

#.. Plotar o histograma dos valores t observados
par(mfrow=c(1, 1))
library(fdth)
plot(fdt(tm,
         k=k),
     type='d',
     xlim=c(-10, 10),
     ylim=c(0, .45),
     main="Histograma e FDP",
     xlab='t',
     ylab="F(t)")

#.. Adicionar a curva teórica dos valores t
# x = t
curve(dt(x,
         df=n-1),
      add=T,
      col="red",
      lwd=3,
      n=1e3)

# Adicionar os limites para a inferência
qtl = qt(p/2, df=n-1, lower.tail=T)
qtu = qt(p/2, df=n-1, lower.tail=F)

segments(qtl, -1, qtl, 0.4, col='red', lwd=1, lty=2)
segments(qtu, -1, qtu, 0.4, col='red', lwd=1, lty=2)

text(-5, 0.30, 'RRHo', col='red')
text(-5, 0.28, paste(100 * p/2, '%', sep=''), cex=.9, col='red')

text( 5, 0.30, 'RRHo', col='red')
text( 5, 0.28, paste(100 * p/2, '%', sep=''), cex=.9, col='red')

text( 0, 0.43, 'RAHo', col='darkgreen')

text( 6, 0.15, 'H0: t=0', pos=4)
text( 6, 0.12, 'H1: t != 0', pos=4)

#.. Adicionar a curva teórica normal reduzida
curve(dnorm(x,
            mean=0,
            sd=1),
      add=T,
      col=gray(.5),
      lwd=2,
      n=1e3)