Processos INAR (1)

Processo INAR (1)

Objetivo

• Utilizar procesos INAR (1), tendo como base três distribuições diferentes no termo aleatório do erro. E em nosso presente estudo as distribuições são Poisson, Binomial Negativa e Geométrica.

• Utilizar o operador thinning binomial, que utiliza da distribuição binomial.

• A forma do processo pode ser vista da seguinte forma:

\[Y_t = \alpha \hspace{0.2cm} \circ \hspace{0.2cm} Y_{t-1} \hspace{0.2cm} + \hspace{0.2cm} \epsilon_t\]

Em que o operador \(\alpha \hspace{0.2cm} \circ = \sum_{i=1}^{Y} N_i\)

• Vamos fazer algumas simulações de séries e testar as estimativas dos parâmetros, através de dois métodos: Mínimos quadrados condicionais e Máxima Verossimilhança Condicional.

• Após a estimativa dos parâmetros vamos simular a série com os parâmetros estimados e comparar com a série verdadeira, e no final vamos utilizar uma série de valores de contagem real (do mundo prático), e submeter aos métodos para observar qual se adequará melhor aos meus dados.

• Observação: Na utilização das distribuições binomial negativa e geométrica utilizaremos ambas parametrizadas pela média, no caso da distribuição binomial negativa, ela será parametrizada por \(\mu\) e por \(\sigma^2\).

Bibliotecas

• De início utilizarei apenas poucas bibliotecas.

• optimx: Utilizada para maximizar a função de log verossimilhança na estimativa dos parâmetros

• DataExplorer: Utilizado para fazer a estatística Descritiva da Série Temporal

• dygraphs: Plotar a série temporal

• plotly: Outra biblioteca para plotar a série temporal

• skimr: Estatística Descritiva da Série Temporal, seja as simuladas ou conjunto de dados reais

• Criaremos aqui em blocos as funções de simular a série temporal, com base nas distribuições, bem como a de estimar usando mínimos quadrados condicionais e máxima verossimilhança condicional.

Simular Série

Code

simular_poisson <- function(n, alpha, lambda){
  x <- numeric(n)
  x[1] <- rpois(1, lambda)
  
  for(t in 2:n){
    thinnig_binomial <- rbinom(1, x[t-1], alpha)
    erro <- rpois(1, lambda)
    x[t] <- thinnig_binomial + erro
  }
  return(x)
}

simular_binomial_negativa <- function(n, alpha, mu, sigma2){
  x <- numeric(n)
  r <- mu^2 / (sigma2 - mu)
  p <- mu / sigma2
  x[1] <- rnbinom(1, size = r, prob = p)
  
  for (t in 2:n) {
    thinning_binomial <- rbinom(1, x[t-1], alpha)
    erro <- rnbinom(1, size = r, prob = p)
    x[t] <- thinning_binomial + erro
  }
  return(x)
}

simular_geometrica <- function(n, alpha, mu){
  x <- numeric(n)
  p <- 1 / (1 + mu)
  x[1] <- rgeom(1, p)
  
  for (t in 2:n) {
    thinning_binomial <- rbinom(1, x[t-1], alpha)
    erro <- rgeom(1, p)
    x[t] <- thinning_binomial + erro
  }
  return(x)
}

Mínimos Quadrados Condicionais

Code

minimos_quadrados <- function(x){
  n <- length(x)
  S1 <- sum(x[-1] * x[-n])                                  # Soma de Y_t * Y_{t-1}
  S2 <- sum(x[-n]^2)                                        # Soma de Y_{t-1}^2
  S3 <- sum(x[-1])                                          # Soma de Y_t
  S4 <- sum(x[-n])                                          # Soma de Y_{t-1}
  
  # Estimado de alpha 
  
  alpha_hat <- (S1 - (S3 * S4) / (n-1)) / (S2 - (S4^2) / (n-1))
  alpha_hat <- max(0, min(alpha_hat, 1)) # Por definição sabemos que alpha pertence ao intervalo [0,1]
  
  lambda_hat <- (S3 - alpha_hat * S4) / (n-1)
  lambda_hat <- max(0, lambda_hat) # Como valores de contagem são positivos, temos que definir esse intervalo
  
  return(c(alpha_estimado = alpha_hat, lambda_estimado = lambda_hat))
}

Máxima Verossimilhança Condicional

Code

logveros_poisson <- function(parametros, x){
  alpha <- parametros[1]
  lambda <- parametros[2]
  n <- length(x)
  log_lik <- 0
  
  for(t in 2:n){
    probabilidade <- 0
    for (k in 0:min(x[t], x[t-1])) {
      probabilidade <- probabilidade + dbinom(k, x[t-1], alpha) * dpois(x[t] - k, lambda)
    }
    log_lik <- log_lik + log(probabilidade)
  }
  return(-log_lik)
}

estimar_max_ver_poisson <- function(x){
  chutes <- c(alpha = 0.4, lambda = mean(x))
  
  resultado <- optimx(
    par = chutes, 
    fn = logveros_poisson, 
    x = x, 
    method = "L-BFGS-B",
    lower = c(0, 0), 
    upper = c(1, Inf)
  )
  return(c(alpha = resultado$alpha, lambda = resultado$lambda))
}

logveros_binomial_negativa <- function(parametros, x){
  alpha <- parametros[1] 
  mu <- parametros[2]
  sigma2 <- parametros[3]
  
  r <- mu^2 / (sigma2 - mu)
  p <- mu / sigma2
  n <- length(x)
  log_lik <- 0
  
  for (t in 2:n) {
    probabilidade <- 0
    for (k in 0:min(x[t], x[t-1])) {
      probabilidade <- probabilidade + dbinom(k, x[t-1], alpha) * dnbinom(x[t] - k , size = r, prob = p)
    }
    log_lik <- log_lik + log(probabilidade)
  }
  return(-log_lik)
}

estimar_max_ver_binomial_negativa <- function(x){
  chutes <- c(alpha = 0.4, mu = mean(x), sigma2 = var(x))
  
  resultado <- optimx(
    par = chutes, 
    fn = logveros_binomial_negativa, 
    x = x, 
    method = "L-BFGS-B", 
    lower = c(0, 0, 0), 
    upper = c(1, Inf, Inf)
  )
  return(c(alpha = resultado$alpha, mu = resultado$mu, sigma = resultado$sigma2))
}

logveros_geometrica <- function(parametros, x){
  alpha <- parametros[1]
  mu <- parametros[2]
  n <- length(x)
  p <- 1 / (1 + mu)
  log_lik <- 0
  
  for(t in 2:n){
    probabilidade <- 0
    for (k in 0:min(x[t], x[t-1])) {
      probabilidade <- probabilidade + dbinom(k, x[t-1], alpha) * dgeom(x[t] - k, p)
    }
    log_lik <- log_lik + log(probabilidade)
  }
  return(-log_lik)
}

estimar_max_ver_geometrica <- function(x){
  chutes <- c(alpha = 0.4, mu = mean(x))
  
  resultado <- optimx(
    par = chutes, 
    fn = logveros_geometrica, 
    x = x, 
    method = "L-BFGS-B", 
    lower = c(0, 0), 
    upper = c(1, Inf)
  )
  return(c(alpha = resultado$alpha, mu = resultado$mu))
}

Simulações

• Vamos simular as séries de tamanhos variados e com as distribuições variadas, sendo variada também as distribuições que usaremos

Observação: Na primeira simulação, vou considerar séries de tamanho 100, com um \(\alpha\) = 0.7 e com \(\lambda\) = 2, no caso da binomial negativa, que temos \(\sigma^2\) = 4.

Code

set.seed(2025)
alpha_verdadeiro <- 0.7
lambda_verdadeiro <- 2
n <- 100
serie_poisson_100 <- simular_poisson(n, alpha = alpha_verdadeiro, lambda = lambda_verdadeiro)
serie_geometrica_100 <- simular_geometrica(n, alpha = alpha_verdadeiro, mu = lambda_verdadeiro)
serie_binomial_negativa_100 <- simular_binomial_negativa(n, alpha = alpha_verdadeiro, mu = lambda_verdadeiro, sigma2 = 4)

Vamos usar os mesmos parâmetros, vamos apenas aumentar o tamanho da série para 150

Code

serie_poisson_150 <- simular_poisson(n = 150, alpha = alpha_verdadeiro, lambda = lambda_verdadeiro)
serie_geometrica_150 <- simular_geometrica(n = 150, alpha = alpha_verdadeiro, mu = lambda_verdadeiro)
serie_binomial_negativa_150 <- simular_binomial_negativa(n = 150, alpha = alpha_verdadeiro, mu = lambda_verdadeiro, sigma2 = 4)

Aumentando para 200

Code

serie_poisson_200 <- simular_poisson(n = 200, alpha = alpha_verdadeiro, lambda = lambda_verdadeiro)
serie_geometrica_200 <- simular_geometrica(n = 200, alpha = alpha_verdadeiro, mu = lambda_verdadeiro)
serie_binomial_negativa_200 <- simular_binomial_negativa(n = 200, alpha = alpha_verdadeiro, mu = lambda_verdadeiro, sigma2 = 4)

Visualizando as séries

Tamanho 100

Code

poisson_100 <- ts(serie_poisson_100)
binomial_negativa_100 <- ts(serie_binomial_negativa_100)
geometrica_100 <- ts(serie_geometrica_100)

series_100 <- cbind(poisson_100, binomial_negativa_100, geometrica_100)
colnames(series_100) <- c("Poisson", "Binomial Negativa", "Geometrica")

dygraph(series_100, main = "Visulização das Séries Temporais de Tamanho 100") %>% 
  dySeries("Poisson", label = "Poisson", drawPoints = TRUE, pointSize = 3) %>% 
  dySeries("Binomial Negativa", label = "Binomial Negativa", drawPoints = TRUE, pointSize = 3) %>% 
  dySeries("Geometrica", label = "Geométrica", drawPoints = TRUE, pointSize = 3 )

Tamanho 150

Code

poisson_150 <- ts(serie_poisson_150)
binomial_negativa_150 <- ts(serie_binomial_negativa_150)
geometrica_150 <- ts(serie_geometrica_150)

series_150 <- cbind(poisson_150, binomial_negativa_150, geometrica_150)
colnames(series_150) <- c("Poisson", "Binomial Negativa", "Geometrica")

dygraph(series_150, main = "Visualização da Série Temporal de tamanho 150") %>% 
  dySeries("Poisson", label = "Poisson", drawPoints = TRUE, pointSize = 3) %>% 
  dySeries("Binomial Negativa", label = "Binomial Negativa", drawPoints = TRUE, pointSize = 3) %>% 
  dySeries("Geometrica", label = "Geométrica", drawPoints = TRUE, pointSize = 3)

Tamanho 200

Code

poisson_200 <- ts(serie_poisson_200)
binomial_negativa_200 <- ts(serie_binomial_negativa_200)
geometrica_200 <- ts(serie_geometrica_200)

series_200 <- cbind(poisson_200, binomial_negativa_200, geometrica_200)
colnames(series_200) <- c("Poisson", "Binomial Negativa", "Geometrica")

dygraph(series_200, main = "Visualização da Série Temporal de tamanho 200") %>% 
  dySeries("Poisson", label = "Poisson", drawPoints = TRUE, pointSize = 3) %>% 
  dySeries("Binomial Negativa", label = "Binomial Negativa", drawPoints = TRUE, pointSize = 3) %>% 
  dySeries("Geometrica", label = "Geométrica", drawPoints = TRUE, pointSize = 3)

Estimando os parâmetros

Vamos estimar os parâmetros e após isso simular a série usando os parâmetros estimados para observar o desempenhos dos métodos.

Code

min_quad_poi_100 <- minimos_quadrados(serie_poisson_100)
min_quad_bin_negativa_100 <- minimos_quadrados(serie_binomial_negativa_100)
min_quad_geometrica_100 <- minimos_quadrados(serie_geometrica_100)


min_quad_poi_150 <- minimos_quadrados(serie_poisson_150)
min_quad_bin_negativa_150 <- minimos_quadrados(serie_binomial_negativa_150)
min_quad_geometrica_150 <- minimos_quadrados(serie_geometrica_150)

min_quad_poi_200 <- minimos_quadrados(serie_poisson_200)
min_quad_bin_negativa_200 <- minimos_quadrados(serie_binomial_negativa_200)
min_quad_geometrica_200 <- minimos_quadrados(serie_geometrica_200)

Code

max_vero_poisson_100 <- estimar_max_ver_poisson(serie_poisson_100)
max_vero_bin_nega_100 <- estimar_max_ver_binomial_negativa(serie_binomial_negativa_100)
max_vero_geom_100 <- estimar_max_ver_geometrica(serie_geometrica_100)

Code

max_vero_poisson_150 <- estimar_max_ver_poisson(serie_poisson_150)
max_vero_bin_nega_150 <- estimar_max_ver_binomial_negativa(serie_binomial_negativa_150)
max_vero_geom_150 <- estimar_max_ver_geometrica(serie_geometrica_150)

Code

max_vero_poisson_200 <- estimar_max_ver_poisson(serie_poisson_200)
max_vero_bin_nega_200 <- estimar_max_ver_binomial_negativa(serie_binomial_negativa_200)
max_vero_geom_200 <- estimar_max_ver_geometrica(serie_geometrica_200)

:::

• Simulação de Monte Carlo

• Fazer a simulação gerando 5.000 séries diferentes

Code

simulacao_monte_carlo <- function(n_simulacoes, n_amostra, alpha_verdadeiro, mu_verdadeiro, sigma_verdadeiro, distribuicao){
  
  variavel <- ifelse(distribuicao %in% c("poisson", "geometrica"), 2, 3)
  resultados_CLS <- matrix(NA, nrow = n_simulacoes, ncol = 2)
  resultados_CML <- matrix(NA, nrow = n_simulacoes, ncol = variavel)
  
  if(distribuicao == "poisson"){
    for (i in 1:n_simulacoes) {
      x <- simular_poisson(n = n_amostra, alpha = alpha_verdadeiro, lambda = mu_verdadeiro) 
    
      resultados_CLS[i, ] <- minimos_quadrados(x)
    
      resultados_CML[i, ] <- estimar_max_ver_poisson(x)
    }
  }
  if(distribuicao == "binomial"){
    for (i in 1:n_simulacoes) {
      x <- simular_binomial_negativa(n = n_amostra, alpha = alpha_verdadeiro, mu = mu_verdadeiro, sigma2 = sigma_verdadeiro)
      resultados_CLS[i, ] <- minimos_quadrados(x)
      
      resultados_CML[i, ] <- estimar_max_ver_binomial_negativa(x)
    }
  }
  if(distribuicao == "geometrica"){
    for (i in 1:n_simulacoes) {
      
      x <- simular_geometrica(n = n_amostra, alpha = alpha_verdadeiro, mu = mu_verdadeiro)
      
      resultados_CLS[i, ] <- minimos_quadrados(x)
      resultados_CML[i, ] <- estimar_max_ver_geometrica(x)
      
    }
  }
  
  media_CLS <- colMeans(resultados_CLS, na.rm = T)
  variancia_CLS <- apply(resultados_CLS, 2, var, na.rm = T)
  
  media_CML <- colMeans(resultados_CML, na.rm = T)
  variancia_CML <- apply(resultados_CML, 2, var, na.rm = T)
  
  return(list(
    media_Minimos_Quadrados_condicionais = media_CLS, 
    variancia_Minimos_Quadrados_condicionais = variancia_CLS, 
    media_Maxima_verossimilhanca_condicional = media_CML, 
    variancia_Maxima_verossimilhanca_condicional = variancia_CML,
    minimos_quadrados_resultado <- resultados_CLS, 
    maxima_verossimilhanca_resultado <- resultados_CML
  ))
}

Simulações

Aqui, inicialmente a primeira simulação teremos os seguintes valores de \(\alpha = 0.6\) e \(\lambda = 4\), e distribuição Geométrica.

Amostra de tamanho 50

Code

resultados_geom_50 <- simulacao_monte_carlo(n_simulacoes = 5000, n_amostra = 50, alpha_verdadeiro = 0.6, mu_verdadeiro = 4, distribuicao = "geometrica")

Code

paste0("Com tamanho de amostra 50 e via mínimos quadrados condicionais, temos o seguinte alpha médio: ", round(resultados_geom_50$media_Minimos_Quadrados_condicionais[1], 4)," para o lambda Médio: ", round(resultados_geom_50$media_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 50 e via mínimos quadrados condicionais, temos o seguinte alpha médio: 0.5451 para o lambda Médio: 4.5356"

Code

paste0("Com tamanho de amostra 50 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: ",   round(resultados_geom_50$media_Maxima_verossimilhanca_condicional[1], 4)," para o lambda Médio: ", round(resultados_geom_50$media_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 50 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: 0.5953 para o lambda Médio: 4.0678"

Code

paste0("Com tamanho de amostra 50 e via Mínimos Quadrados Condicionais, temos a variância do alpha: ", 
       round(resultados_geom_50$variancia_Minimos_Quadrados_condicionais[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_50$variancia_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 50 e via Mínimos Quadrados Condicionais, temos a variância do alpha: 0.0136 e a variância para o lambda: 1.824"

Code

paste0("Com tamanho de amostra 50 e via Máxima Verossimilhança Condicional, temos a variância do alpha: ", 
       round(resultados_geom_50$variancia_Maxima_verossimilhanca_condicional[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_50$variancia_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 50 e via Máxima Verossimilhança Condicional, temos a variância do alpha: 0.0018 e a variância para o lambda: 0.4701"

Code

banco_geom_50_min_quad <- as.data.frame(resultados_geom_50[[5]])
banco_geom_50_max_vero <- as.data.frame(resultados_geom_50[[6]])

• Visualizando os resultados de Maneira Gráfica

• Através de um histograma, podemos analisar quais foram os resultados da estimação dos parâmetros.

Code

grafico_50_min_quad_lambda <- ggplot(banco_geom_50_min_quad) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 fill = "gray", alpha = 0.5) +
 labs(x = "lambda", y = "Quantidade", title = "Histograma da Estimação de lambda com tamanho de amostra 50", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_50_min_quad_alpha <- ggplot(banco_geom_50_min_quad) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "Alpha", y = "Quantidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 50", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_50_max_vero_lambda <- ggplot(banco_geom_50_max_vero) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "lambda", y = "Quantidade", title = "Histograma da Estimação de lambda com tamanho de amostra 50", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_50_max_vero_alpha <- ggplot(banco_geom_50_max_vero) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "Alpha", y = "Quantidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 50", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

gridExtra::grid.arrange(grafico_50_min_quad_alpha, grafico_50_max_vero_alpha, grafico_50_min_quad_lambda, grafico_50_max_vero_lambda, nrow = 2, ncol = 2)

Amostra de tamanho 100

Code

resultados_geom_100 <- simulacao_monte_carlo(n_simulacoes = 5000, n_amostra = 100, alpha_verdadeiro = 0.6, mu_verdadeiro = 4, distribuicao = "geometrica")

Code

paste0("Com tamanho de amostra 100 e via mínimos quadrados condicionais, temos o seguinte alpha médio: ", round(resultados_geom_100$media_Minimos_Quadrados_condicionais[1], 4)," para o lambda Médio: ", round(resultados_geom_100$media_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 100 e via mínimos quadrados condicionais, temos o seguinte alpha médio: 0.5712 para o lambda Médio: 4.2801"

Code

paste0("Com tamanho de amostra 100 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: ",   round(resultados_geom_100$media_Maxima_verossimilhanca_condicional[1], 4)," para o lambda Médio: ", round(resultados_geom_100$media_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 100 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: 0.5928 para o lambda Médio: 4.0768"

Code

paste0("Com tamanho de amostra 100 e via Mínimos Quadrados Condicionais, temos a variância do alpha: ", 
       round(resultados_geom_100$variancia_Minimos_Quadrados_condicionais[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_100$variancia_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 100 e via Mínimos Quadrados Condicionais, temos a variância do alpha: 0.0067 e a variância para o lambda: 0.8526"

Code

paste0("Com tamanho de amostra 100 e via Máxima Verossimilhança Condicional, temos a variância do alpha: ", 
       round(resultados_geom_100$variancia_Maxima_verossimilhanca_condicional[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_100$variancia_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 100 e via Máxima Verossimilhança Condicional, temos a variância do alpha: 9e-04 e a variância para o lambda: 0.2573"

Code

banco_geom_100_min_quad <- as.data.frame(resultados_geom_100[[5]])
banco_geom_100_max_vero <- as.data.frame(resultados_geom_100[[6]])

• Visualizando os resultados de Maneira Gráfica

• Através de um histograma, podemos analisar quais foram os resultados da estimação dos parâmetros.

Code

grafico_100_min_quad_lambda <- ggplot(banco_geom_100_min_quad) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "lambda", y = "Densidade", title = "Histograma da Estimação de lambda com tamanho de amostra 100", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_100_min_quad_alpha <- ggplot(banco_geom_100_min_quad) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density(
 alpha = 0.5, fill = "gray") +
 labs(x = "Alpha", y = "Densidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 100", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_100_max_vero_lambda <- ggplot(banco_geom_100_max_vero) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "lambda", y = "Densidade", title = "Histograma da Estimação de lambda com tamanho de amostra 100", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_100_max_vero_alpha <- ggplot(banco_geom_100_max_vero) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density(
 alpha = 0.5, fill = "gray") +
 labs(x = "Alpha", y = "Densidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 100", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

gridExtra::grid.arrange(grafico_100_min_quad_alpha, grafico_100_max_vero_alpha, grafico_100_min_quad_lambda, grafico_100_max_vero_lambda, nrow = 2, ncol = 2)

Amostra de tamanho 300

Code

resultados_geom_300 <- simulacao_monte_carlo(n_simulacoes = 5000, n_amostra = 300, alpha_verdadeiro = 0.6, mu_verdadeiro = 4, distribuicao = "geometrica")

Code

paste0("Com tamanho de amostra 300 e via mínimos quadrados condicionais, temos o seguinte alpha médio: ", round(resultados_geom_300$media_Minimos_Quadrados_condicionais[1], 4)," para o lambda Médio: ", round(resultados_geom_300$media_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 300 e via mínimos quadrados condicionais, temos o seguinte alpha médio: 0.5909 para o lambda Médio: 4.0918"

Code

paste0("Com tamanho de amostra 300 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: ",   round(resultados_geom_300$media_Maxima_verossimilhanca_condicional[1], 4)," para o lambda Médio: ", round(resultados_geom_300$media_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 300 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: 0.6114 para o lambda Médio: 3.9755"

Code

paste0("Com tamanho de amostra 300 e via Mínimos Quadrados Condicionais, temos a variância do alpha: ", 
       round(resultados_geom_300$variancia_Minimos_Quadrados_condicionais[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_300$variancia_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 300 e via Mínimos Quadrados Condicionais, temos a variância do alpha: 0.0022 e a variância para o lambda: 0.2849"

Code

paste0("Com tamanho de amostra 300 e via Máxima Verossimilhança Condicional, temos a variância do alpha: ", 
       round(resultados_geom_300$variancia_Maxima_verossimilhanca_condicional[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_300$variancia_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 300 e via Máxima Verossimilhança Condicional, temos a variância do alpha: 4e-04 e a variância para o lambda: 0.0915"

Code

banco_geom_300_min_quad <- as.data.frame(resultados_geom_300[[5]])
banco_geom_300_max_vero <- as.data.frame(resultados_geom_300[[6]])

• Visualizando os resultados de Maneira Gráfica

• Através de um histograma, podemos analisar quais foram os resultados da estimação dos parâmetros.

Code

grafico_300_min_quad_lambda <- ggplot(banco_geom_300_min_quad) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "lambda", y = "Densidade", title = "Histograma da Estimação de lambda com tamanho de amostra 300", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_300_min_quad_alpha <- ggplot(banco_geom_300_min_quad) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "Alpha", y = "Densidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 300", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_300_max_vero_lambda <- ggplot(banco_geom_300_max_vero) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "lambda", y = "Densidade", title = "Histograma da Estimação de lambda com tamanho de amostra 300", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_300_max_vero_alpha <- ggplot(banco_geom_300_max_vero) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density(
 alpha = 0.5, fill = "#112446") +
 labs(x = "Alpha", y = "Densidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 300", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

gridExtra::grid.arrange(grafico_300_min_quad_alpha, grafico_300_max_vero_alpha, grafico_300_min_quad_lambda, grafico_300_max_vero_lambda, nrow = 2, ncol = 2)

Amostra de tamanho 500

Code

resultados_geom_500 <- simulacao_monte_carlo(n_simulacoes = 5000, n_amostra = 500, alpha_verdadeiro = 0.6, mu_verdadeiro = 4, distribuicao = "geometrica")

Code

paste0("Com tamanho de amostra 500 e via mínimos quadrados condicionais, temos o seguinte alpha médio: ", round(resultados_geom_500$media_Minimos_Quadrados_condicionais[1], 4)," para o lambda Médio: ", round(resultados_geom_500$media_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 500 e via mínimos quadrados condicionais, temos o seguinte alpha médio: 0.5941 para o lambda Médio: 4.0603"

Code

paste0("Com tamanho de amostra 500 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: ",   round(resultados_geom_500$media_Maxima_verossimilhanca_condicional[1], 4)," para o lambda Médio: ", round(resultados_geom_500$media_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 500 e via Máxima Verossimilhança Condicional, temos o seguinte alpha médio: 0.6126 para o lambda Médio: 3.9803"

Code

paste0("Com tamanho de amostra 500 e via Mínimos Quadrados Condicionais, temos a variância do alpha: ", 
       round(resultados_geom_500$variancia_Minimos_Quadrados_condicionais[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_500$variancia_Minimos_Quadrados_condicionais[2], 4))

[1] "Com tamanho de amostra 500 e via Mínimos Quadrados Condicionais, temos a variância do alpha: 0.0014 e a variância para o lambda: 0.1753"

Code

paste0("Com tamanho de amostra 500 e via Máxima Verossimilhança Condicional, temos a variância do alpha: ", 
       round(resultados_geom_500$variancia_Maxima_verossimilhanca_condicional[1], 4), " e a variância para o lambda: ", 
       round(resultados_geom_500$variancia_Maxima_verossimilhanca_condicional[2], 4))

[1] "Com tamanho de amostra 500 e via Máxima Verossimilhança Condicional, temos a variância do alpha: 2e-04 e a variância para o lambda: 0.057"

Code

banco_geom_500_min_quad <- as.data.frame(resultados_geom_500[[5]])
banco_geom_500_max_vero <- as.data.frame(resultados_geom_500[[6]])

• Visualizando os resultados de Maneira Gráfica

• Através de um histograma, podemos analisar quais foram os resultados da estimação dos parâmetros.

Code

grafico_500_min_quad_lambda <- ggplot(banco_geom_500_min_quad) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "lambda", y = "Densidade", title = "Histograma da Estimação de lambda com tamanho de amostra 500", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_500_min_quad_alpha <- ggplot(banco_geom_500_min_quad) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "Alpha", y = "Densidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 500", 
      subtitle = "Via Mínimos Quadrados Condicionais") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_500_max_vero_lambda <- ggplot(banco_geom_500_max_vero) +
 aes(x = V2) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "lambda", y = "Densidade", title = "Histograma da Estimação de lambda com tamanho de amostra 500", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

grafico_500_max_vero_alpha <- ggplot(banco_geom_500_max_vero) +
 aes(x = V1) +
 geom_histogram(aes(y = ..density..), bins = 60L, fill = "#46B47A", color = "black") +
 geom_density( 
 alpha = 0.5, fill = "gray") +
 labs(x = "Alpha", y = "Densidade", title = "Histograma da Estimação de Alpha com tamanho de amostra 500", 
      subtitle = "Via Máxima Verossimilhança Condicional") +
 theme_classic() +
 theme(plot.title = element_text(size = 12L, hjust = 0.5))

Code

gridExtra::grid.arrange(grafico_500_min_quad_alpha, grafico_500_max_vero_alpha, grafico_500_min_quad_lambda, grafico_500_max_vero_lambda, nrow = 2, ncol = 2)