#===============================================================================
# Título: Inferência Estatística - Entender origem da distribuição Qui-Quadrado
# Autor : José Cláudio Faria/UESC/DCET
# Data  : 21/08/2019 12:57:10
# Versão: v6
#===============================================================================
# Objetivos
#===============================================================================
# a) Gerar uma população normal;
# b) Amostrar repetidamente a população gerada;
# c) Em cada amostra: elevar ao quadrado cada elemento e realizar o somatório
#    armazenando o valor num vetor (qq);
# d) Plotar o histograma dos valores (qq), que deve ser semelhante a função X2,
#    para o tamanho da amostra, 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 associada ao
#    tamanho da amostra (n-1).
#===============================================================================
# Link com boas informações: # http://www.portalaction.com.br/content/63-distribui%C3%A7%C3%A3o-qui-quadrado
#===============================================================================

#! CASO 1: A PARTIR DA NORMAL PADRÃO (usando a média populacional no cálculo do desvio)
#!. Dados
n  <- 5   # Tamanho da amostra
k  <- 100 # Número de classes do histograma
nr <- 1e4 # Número de repetições do processo de moastragem

# Amostragem
m <- matrix(rnorm(nr*n),
            nc=n)
str(m)
head(m)
# Cada linha de m pode ser considerada como uma amostra de tamanho 'n' de POP.
# POP ~ N (0, 1)

qq <- apply(m,
            1,
            function(x) sum(x^2))
# O vetor 'qq' deve apresentar distribuição Qui-Quadrado com n gl!

library(fdth)
par(mfrow=c(2, 2))
plot(fdt(qq,
         k=k),
     typ='d',
     xlab='',
     main='Origem: N(0, 1)')

# df=n-1: não explica muito bem!
curve(dchisq(x,
             df=n-1),
      lw=3,
      col=gray(.5),
      add=T,
      n=1e3)

# df=n: perfeito!
curve(dchisq(x,
             df=n),
      lw=3,
      col='darkblue',
      add=T,
      n=1e3)


#! CASO 2: A PARTIR DA NORMAL PADRÃO (usando a estimativa da média no cálculo do desvio)
#!. Dados
n  <- 5   # Tamanho da amostra
k  <- 100 # Número de classes do histograma
nr <- 1e4 # Número de repetições do processo de moastragem

# Amostragem
m <- matrix(rnorm(nr*n),
            nc=n)
str(m)
head(m)
# Cada linha de m pode ser considerada como uma amostra de tamanho 'n' de POP.
# POP ~ N (0, 1)

qq <- apply(m,
            1,
            function(x) sum((x - mean(x))^2))
# O vetor 'qq' deve apresentar distribuição Qui-Quadrado com n-1 gl!

plot(fdt(qq,
         k=k),
     typ='d',
     xlab='',
     main='Origem: N(0, 1)')

curve(dchisq(x,
             df=n-1),
      lw=3,
      col='darkblue',
      add=T,
      n=1e3)


#! CASO 3: GENERALIZA PARA QUALQUER NORMAL (usando a estimativa da média no cálculo do desvio)
#!. Dados
m  <- 10  # Média de POP
sd <-  2  # Desvio padrão de POP
n  <-  5  # Tamanho da amostra
k  <- 100 # Número de classes do histograma
nr <- 1e4 # Número de repetições do processo de moastragem

# Amostragem
m <- matrix(rnorm(nr*n,
                  m,
                  sd),
            nc=n)
str(m)
head(m)
# Cada linha de m pode ser considerada como uma amostra de tamanho 'n' de POP.
# POP ~ N (0, 1)

qq <- apply(m,
            1,
            function(x) sum((x - mean(x))^2)/sd^2)  # Faz a padronização!
# O vetor 'qq' deve apresentar distribuição Qui-Quadrado com n-1 gl!

plot(fdt(qq,
         k=k),
     typ='d',
     xlab='',
     main='Origem: N(10, 2)')

curve(dchisq(x,
             df=n-1),
      lw=3,
      col='darkblue',
      add=T,
      n=1e3)


#! CASO 4: GENERALIZA PARA QUALQUER NORMAL (usando a estimativa da média no cálculo do desvio)
#!. Dados
m  <- 20  # Média de POP
sd <-  4  # Desvio padrão de POP
n  <-  5  # Tamanho da amostra
k  <- 100 # Número de classes do histograma
nr <- 1e4 # Número de repetições do processo de moastragem

# Amostragem
m <- matrix(rnorm(nr*n,
                  m,
                  sd),
            nc=n)
str(m)
head(m)
# Cada linha de m pode ser considerada como uma amostra de tamanho 'n' de POP.
# POP ~ N (0, 1)

qq <- apply(m,
            1,
            function(x) (n-1)*var(x)/sd^2)
# O vetor 'qq' deve apresentar distribuição Qui-Quadrado com n-1 gl!

plot(fdt(qq,
         k=k),
     typ='d',
     main='Origem: N(20, 4)')

curve(dchisq(x,
             df=n-1),
      lw=3,
      col='darkblue',
      add=T,
      n=1e3)