#===============================================================================
# Name   : Comparação das estimativas da medidas estatísticas de dispersão:
#          quadrado médio dos desvios (qmd) e variância (s2)
# Author : J.C.Faria
# Date   : 2024-11-14 18:26:03
# Aim    : Ensino e percepção dos fundamentos das medida estatísticas
#          de dispersão
# Mail   : <<<joseclaudio.faria@gmail.com>>>
#===============================================================================

#.. ___Ini opções___
#... Tamanho da amostra
n <- 5      # Usar >= 3

#... Características da população
Mu     <- 10
Sigma2 <- 2

#... Número de amostras
k <- 200

#...limites do eixo Y
y_min <- 0
y_max <- 2*Sigma2
#.. ___Fim opções___

#. Simulação
#.. Funções
# Estimativa dos quadrado médio dos desvios - qmd
f_qmd <- function(x)
{
  sum((x - mean(x))^2)/length(x)
}

# Estimativa da variância
f_s2 <- function(x)
{
  sum((x - mean(x))^2)/(length(x) - 1)
}

#.. Gerando população
pop <- rnorm(k*n,
             Mu,
             sqrt(Sigma2))
                          
#.. Estimativas da variância: f_qmd
qmd <- apply(matrix(pop,
                    ncol=n),
             1,
             f_qmd)

#.. Estimativas da variância: f_s2
s2 <- apply(matrix(pop,
                   ncol=n),
            1,
            f_s2)
               

# Esperanças
E_qmd  <- mean(qmd)
E_s2   <- mean(s2)


#. Plot
layout(matrix(1:2,
              ncol=2),
       widths=c(1,
                3),
       heights=c(3.5,
                 3.5),
       respect=FALSE)

boxplot(pop,
        main= 'POP',
        col='gold')

plot(1:4,
     xaxt='n',
     yaxt='n',
     ylab='Parâmetro e estimativas',
     ylim=c(y_min, 
            y_max),
     type='n',
     xlab='',
     main=paste('Tamanho da amostra:',
                ' n=',   
                 n,
                 sep=''))

axis(s=2,
     at=Sigma2,
     labels=expression('QMD' == sigma^2))

axis(s=1,
     at=c(2.5, 3.5),
     labels=c(k, k),
     font=3)

#... A variância informada
abline(h=Sigma2,
       lwd=3)

#... Esperança da estimativa da variância: f_s2
points(x=rep(2.5, k),
       y=s2,
       pch='.',
       col=2,
       cex=2)                 
abline(h=E_s2,
       col=2,
       lwd=2)
text(x=2.5,
     y=E_s2,
     labels=expression(E(s^2)),
     pos=3,
     offset=.2,
     col=2,
     cex=1.5)

#... Esperança da estimativa: f_dqm
points(x=rep(3.5, k),
       y=qmd,
       pch='.',
       col=4,
       cex=2)                 
abline(h=E_qmd,
       col=4,
       lwd=2)
text(x=3.5,
     y=E_qmd,
     labels='E(qmd)',
     pos=3,
     offset=.2,
     col=4,
     cex=1.5)

legend('topleft',
       legend=c(expression(sigma^2 == sum(Di^2, i==1, n)/N),
                expression(s^2 == sum(di^2, i==1, n)/(n-1)),
                expression(qmd == sum(di^2, i==1, n)/n)),
       lwd=c(3, 2, 2),
       col=c(1, 2, 4))

legend('bottomleft',
       legend=c(expression(bold(D[i] == (y[i] - mu))),
                expression(bold(d[i] == (y[i] - m)))))

#.. Identificação no gráfico
# função: f_s2
text(x=2.5,
     y=0.1,
     labels='f_s2',
     pos=1,
     offset=.1,
     col=2,
     cex=1.5)

# função: qmd
text(x=3.5,
     y=0.1,
     labels='f_qmd',
     pos=1,
     offset=.1,
     col=4,
     cex=1.5)