Introdução

A Mega-Sena é a maior modalidade lotérica atualmente no Brasil. Lançada em 1996, a loteria é organizada pela Caixa Econômica Federal (CAIXA). Para apostar na Mega-Sena, o jogador escolhe de 6 a 15 números entre as dezenas numeradas de 01 a 60.

Os sorteios das seis dezenas premiadas são realizados duas vezes por semana. Caso acerte os 6 números sorteados, o jogador leva o prêmio correspondente à Sena (faixa 1). Se acertar 5 números, ganha o prêmio da Quina (faixa 2). Por outro lado, se acertar 4 números, leva o prêmio da Quadra (faixa 3).

Essas informações estão no site de Loterias da CAIXA, onde se pode ver, acessando o link correspondente à Mega-Sena, a seguinte imagem de chamada à aposta:

Chamada à aposta na Mega-Sena

Por outro lado, o próprio site fornece dados valiosos para que o jogador possa tomar suas decisões de aposta. Nesse diapasão, faremos uma análise dos dados disponíveis para a Mega-Sena em busca de responder às seguintes perguntas:

  1. Existe alguma estratégia de jogo mais favorável ao apostador?
  2. É possível determinar a relação custo x retorno das apostas?

Noções elementares sobre probabilidades e loterias

Com o intuito de responder às perguntas levantadas, precisamos compreender algumas premissas, que são características dos jogos de loteria ou consequências da Teoria das Probabilidades :

Em adição às premissas supracitadas, faz-se necessário entender como são calculadas as probabilidades de ganhar em uma loteria tal qual a Mega-Sena.

Um jogo de loteria nada mais é do que a aposta numa combinação de números, com a esperança de que a combinação escolhida seja sorteada para ganhar o prêmio. Esse sorteio favorável possui uma probabilidade de ocorrência, que pode ser mensurada usando os conceitos de Combinação e Probabilidade Frequentista.

Uma combinação é essencialmente a contagem do número de escolhas dentro de um conjunto de possibilidades, independentemente da ordenação dos elementos escolhidos. Na Mega-Sena, uma aposta simples consiste em selecionar exatos 6 números em um conjunto de 60 possíveis escolhas (as dezenas de 01 a 60).

Em suma, o processo de contar quantas escolhas de apostas simples são possíveis na Mega-Sena consiste em enumerar de quantas maneiras o jogador pode preencher 6 espaços vazios em seu cartão com os 60 números disponíveis. Para o primeiro espaço, tem 60 possibilidades, para o segundo 59, para o terceiro 58 e assim por diante. O produto das possibilidades contabiliza o total de escolhas “embaralhadas” dos números. Por causa disso, não se pode esquecer de dividir o total obtido pela quantidade de “embaralhadas” em cada grupo de 6 números escolhidos, já que a permutação na ordem deles leva à contagem de apostas repetidas. Assim, o total de possíveis apostas simples na Mega-Sena é dado matematicamente pela expressão abaixo:

\[ C_6^{60} = \frac{60!}{6!(60-6)!} = \frac{60 \times 59 \times 58 \times 57 \times 56 \times 55}{6 \times 5 \times 4 \times 3 \times 2 \times 1} = 50.063.860 \]

Além da notação de combinação introduzida acima, é necessário entender também a definição clássica de probabilidade. Quando queremos saber a chance de ocorrer um evento, intuitivamente avaliamos os casos que são favoráveis em relação ao total de casos possíveis. Em resumo, essa relação, mais especificamente, essa razão entre os casos favoráveis e os possíveis constitui a chamada definição frequentista de probabilidade.

Ela pode ser usada para verificar, por exemplo, que a probabilidade de acertar a Sena ao jogar um cartão com 8 números é a razão entre o total de casos favoráveis (nesse contexto as combinações de 6 números que estão entre os 8 selecionados) e o total de casos possíveis (todas as combinações de 6 números entre os 60 disponíveis), isto é:

\[ \frac{C_6^{8}}{C_6^{60}} = \frac{28}{50.063.860} = \frac{1}{1.787.995} \]

Unindo esses conceitos a outros de análise combinatória, é possível demonstrar resultados curiosos como:

Tudo isso e mais um pouco está nas entrelinhas dos dados sobre a Mega-Sena. São eles que iremos explorar nesta análise, a fim de determinar se existe uma estratégia vencedora na Mega-Sena que atenda às premissas, com uma relação custo x retorno efetivamente calculável.

Com tal finalidade, usaremos dados das seguintes fontes: Página da Mega-Sena no Site da CAIXA, para obter os dados principais, e Plataforma SIDRA do IBGE, para obter dados da série histórica do Índice Nacional de Preços ao Consumidor (INPC). Ademais, aplicaremos ferramentas de análise de dados da linguagem de programação R, através do Software R versão 4.0.5, na IDE RStudio versão 1.4.1103. Começamos apresentando os pacotes do repositório CRAN do R a serem utilizados, os quais, após instalados pelo comando “install.packages(package_name)”, podem ser carregados da forma a seguir:

# Carregamento dos pacotes a serem utilizados

library(rvest) # Leitura de dados no padrão html
library(sidrar) # Leitura de dados da plataforma SIDRA do IBGE
library(dplyr) # Manipulação de bancos de dados
library(tidyr) # Transformação de dados
library(readr) # Recursos de leitura e formatação de dados
library(stringr) # Manipulação de strings
library(ggplot2) # Recursos gráficos
library(lpSolve) # Otimização linear e inteira

A leitura dos dados pelo R diretamente das fontes indicadas pode ser feita por meio dos seguintes comandos: 

# Leitura dos dados

# Resultados da Mega Sena no site da CAIXA
url_link <-"http://loterias.CAIXA.gov.br/wps/portal/loterias/landing/megasena/!ut/p/a1/04_Sj9CPykssy0xPLMnMz0vMAfGjzOLNDH0MPAzcDbwMPI0sDBxNXAOMwrzCjA0sjIEKIoEKnN0dPUzMfQwMDEwsjAw8XZw8XMwtfQ0MPM2I02-AAzgaENIfrh-FqsQ9wNnUwNHfxcnSwBgIDUyhCvA5EawAjxsKckMjDDI9FQE-F4ca/dl5/d5/L2dBISEvZ0FBIS9nQSEh/pw/Z7_HGK818G0K8DBC0QPVN93KQ10G1/res/id=historicoHTML/c=cacheLevelPage/=/"
resultados <- url_link %>%
  read_html() %>%
  html_nodes(xpath='/html/body/table/tbody')
resultados <- html_table(resultados, header = F, dec = ",")

# Série da variação mensal do INPC desde janeiro/1994 obtida da plataforma SIDRA do IBGE via parâmetros para API
inpc <- get_sidra(api="/t/1736/n1/all/v/44/p/last%20331/d/v44%202")
## All others arguments are desconsidered when 'api' is informed
# Limpeza e tratamento dos dados

resultados <- data.frame(resultados[[1]][,1:21])
colnames(resultados)<-c("Concurso","Local","Data_sorteio","Coluna1","Coluna2","Coluna3","Coluna4","Coluna5","Coluna6",
                                                "Ganhadores_faixa1","Ganhadores_faixa2","Ganhadores_faixa3","Rateio_faixa1","Rateio_faixa2","Rateio_faixa3",
                                                "Cidade","Valor_arrecadado","Estimativa_proximoconcurso","Valoracumulado_proximoconcurso",
                                                "Acumulado","Sorteio_especial")

resultados <- resultados %>% drop_na()

resultados <- resultados %>% mutate(Concurso=as.integer(Concurso),
                      Data_sorteio=as.Date(Data_sorteio,"%d/%m/%Y"),
                      Periodo_referencial=format(Data_sorteio,"%Y-%m"),                    
                      Ano_sorteio=as.integer(format(Data_sorteio,"%Y")),
                      Coluna1=Coluna1,
                      Coluna2=Coluna2,
                      Coluna3=Coluna3,
                      Coluna4=Coluna4,
                      Coluna5=Coluna5,
                      Coluna6=Coluna6,
                      Ganhadores_faixa1=Ganhadores_faixa1,
                      Ganhadores_faixa2=Ganhadores_faixa2,
                      Ganhadores_faixa3=Ganhadores_faixa3,
                      Cidade=Cidade)

resultados <- resultados %>% mutate_at(vars(Rateio_faixa1,Rateio_faixa2,Rateio_faixa3), parse_number, locale = locale(decimal_mark = ","))

resultados <- resultados %>% dplyr::select(1,3,22:23,4:16)

resultados <- resultados %>% filter(Data_sorteio<=as.Date("2021-08-31")) %>% mutate_at(vars(Estado=Cidade), str_sub, -2) %>% dplyr::select(!c(Cidade))

inpc <- inpc %>% mutate(INPC_per100_mensal=Valor) %>% mutate_at(vars(Ref_ano=`Mês (Código)`),str_sub,1,4) %>%
  mutate_at(vars(Ref_mês=`Mês (Código)`),str_sub,-2)
  
inpc <- inpc %>% filter(Ref_ano!="1994"&Ref_ano!="1995") %>% dplyr::select(11:14)

inpc <- inpc %>% mutate(Periodo_referencial=paste0(Ref_ano,"-",Ref_mês)) %>%
  dplyr::select(Periodo_referencial,INPC_per100_mensal)

Após alguns ajustes nos nomes das variáveis e transformações de dados, com os recursos do pacote “tidyr”, obtemos uma melhor visualização das variáveis presentes nos dados e seus respectivos tipos e primeiros valores, por meio da função “glimpse(dataset)” do pacote “dplyr”.

glimpse(resultados)
## Rows: 2,404
## Columns: 17
## $ Concurso            <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,~
## $ Data_sorteio        <date> 1996-03-11, 1996-03-18, 1996-03-25, 1996-04-01, 1~
## $ Periodo_referencial <chr> "1996-03", "1996-03", "1996-03", "1996-04", "1996-~
## $ Ano_sorteio         <int> 1996, 1996, 1996, 1996, 1996, 1996, 1996, 1996, 19~
## $ Coluna1             <int> 4, 9, 10, 1, 1, 7, 3, 4, 8, 4, 15, 4, 18, 2, 12, 2~
## $ Coluna2             <int> 5, 37, 11, 5, 2, 13, 5, 17, 43, 18, 25, 16, 32, 16~
## $ Coluna3             <int> 30, 39, 29, 6, 6, 19, 20, 37, 54, 21, 37, 19, 47, ~
## $ Coluna4             <int> 33, 41, 30, 27, 16, 22, 21, 38, 55, 25, 38, 20, 50~
## $ Coluna5             <int> 41, 43, 36, 42, 19, 40, 38, 47, 56, 38, 58, 27, 54~
## $ Coluna6             <int> 52, 49, 47, 59, 46, 47, 56, 53, 60, 57, 59, 43, 56~
## $ Ganhadores_faixa1   <int> 0, 1, 2, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1,~
## $ Ganhadores_faixa2   <int> 17, 65, 62, 39, 98, 109, 100, 60, 17, 251, 148, 10~
## $ Ganhadores_faixa3   <int> 2016, 4488, 4261, 3311, 5399, 7147, 5736, 5262, 21~
## $ Rateio_faixa1       <dbl> 0.0, 2307162.2, 391192.5, 0.0, 0.0, 0.0, 0.0, 0.0,~
## $ Rateio_faixa2       <dbl> 39158.92, 14424.02, 10515.93, 15322.24, 5318.10, 7~
## $ Rateio_faixa3       <dbl> 330.21, 208.91, 153.01, 180.48, 96.53, 110.03, 152~
## $ Estado              <chr> "", "PR", "RN", "", "", "", "", "", "", "", "RJ", ~


No caso específico do índice de inflação INPC, criamos duas novas colunas a partir da variação percentual mensal do INPC, uma correspondente ao número-índice do mês e outra, ao índice acumulado até a data atual (última divulgação em 10/08/2021, índice referente ao mês de julho).

# Ajustes para obter mais informações sobre o INPC

inpc <- inpc %>% 
  mutate(Num_ind_mensal = 1+INPC_per100_mensal/100) %>%
  mutate(INPC_acum_atual = prod(Num_ind_mensal)/cumprod(Num_ind_mensal)*Num_ind_mensal)

inpc <- rbind(inpc,c("2021-08",0,1,1))

inpc <- inpc %>% mutate(INPC_per100_mensal = round(as.numeric(INPC_per100_mensal),2),
                Num_ind_mensal = round(as.numeric(Num_ind_mensal),4),
                INPC_acum_atual = round(as.numeric(INPC_acum_atual),2))

glimpse(inpc) #Última coluna permite atualizar todos os prêmios passados ao valor presente
## Rows: 311
## Columns: 4
## $ Periodo_referencial <chr> "1996-01", "1996-02", "1996-03", "1996-04", "1996-~
## $ INPC_per100_mensal  <dbl> 1.46, 0.71, 0.29, 0.93, 1.28, 1.33, 1.20, 0.50, 0.~
## $ Num_ind_mensal      <dbl> 1.0146, 1.0071, 1.0029, 1.0093, 1.0128, 1.0133, 1.~
## $ INPC_acum_atual     <dbl> 4.96, 4.89, 4.86, 4.84, 4.80, 4.74, 4.68, 4.62, 4.~

Depois de uma análise exploratória dos dados, usando medidas resumo como média, mediana, quartis, etc., os recursos gráficos, especialmente do pacote “ggplot2”, permitem construir visualizações interessantes. 

barplot(table(c(resultados$Coluna1,resultados$Coluna2,resultados$Coluna3,resultados$Coluna4,resultados$Coluna5,resultados$Coluna6)),width=0.8,col="yellow",main="Frequência dos números sorteados",xlab="Números",ylab="Frequência",xlim=c(1,60),ylim=c(0,300),beside=T)


A distribuição das frequências dos números sorteados não parece contradizer a primeira premissa (sorteios aleatórios), embora o nível de aleatoriedade possa ser submetido a testes que não serão abordados aqui. 

resultados %>% dplyr::select(Ganhadores_faixa1,Estado) %>% filter(Estado!=""&Estado!="--") %>% group_by(Estado) %>%
summarise(Ganhadores_Sena_estado=n()) %>% mutate(Prop_ganhadores_Sena_estado=round(Ganhadores_Sena_estado/sum(Ganhadores_Sena_estado)*100,2)) %>%
  ggplot(aes(Estado,Prop_ganhadores_Sena_estado)) + ggtitle("Localização dos ganhadores") +
  geom_bar(stat="identity",aes(color=Estado),show.legend=FALSE) + geom_label(aes(label = Prop_ganhadores_Sena_estado),vjust = -0.2,angle=90) + 
  ylab("Porcentagem de acertadores da Sena") + theme(plot.title = element_text(hjust = 0.5))


Ao comparar a visualização com dados de Distribuição da População Brasileira, observa-se uma pequena divergência em relação aos estados de origem dos ganhadores da Sena, com peso maior destes nos estados do Centro-Sul, o que indicaria uma tendência maior ao consumo da loteria na referida macrorregião econômica.

Curioso notar que o estado do Amapá (AP) - 2º menor em população no Brasil - é o único que supostamente nunca teve um ganhador do maior prêmio da Mega-Sena. Ressalte-se que eventuais felizardos amapaenses poderiam estar entre os que ganharam por canal eletrônico, representados por XX no gráfico, que não tiveram o estado de origem divulgado. 

resultados %>% dplyr::select(starts_with("Rateio") | ends_with("sorteio")) %>% filter(Rateio_faixa1>0) %>% group_by(Ano_sorteio) %>%
  ggplot(aes(as.factor(Ano_sorteio),Rateio_faixa1/1000000)) + ggtitle("Premiação da Sena ao longo dos anos") +
  geom_boxplot(coef=3) + geom_jitter(width = 0.1, alpha = 0.2) +
  scale_y_continuous(trans="log10") + xlab("Ano do sorteio") +
  ylab("Prêmio (em milhões de R$)") + theme(plot.title = element_text(hjust = 0.5), axis.text.x = element_text(angle=90, hjust=1))


O crescimento dos prêmios, ainda não corrigidos pelo índice de inflação INPC, pode incitar a tentação de querer apostar mais na Mega-Sena hoje, porque ao longo do tempo os quantis das premiações aumentaram de uma maneira geral, apesar de existirem outliers muito baixos (pontos pretos) nos prêmios de 2021.

No entanto, há que se ponderar também a variação do valor do dinheiro no tempo e a própria alteração no custo das apostas (não presente nas fontes de dados consultadas). Essa ausência de dados pode ser mitigada pela atualização dos valores dos prêmios ao presente, com o intuito de que eles possam ser comparados com o custo atual de cada cartão.

quickplot(seq(1996,2021,(2021-1996)/(length(inpc$Num_ind_mensal)-1)),inpc$Num_ind_mensal,xlab="Ano",ylab="INPC mensal",main="Variação do INPC",geom="line", colour = I("blue"))

O gráfico da variação mensal do INPC indica o efeito inflacionário ao longo do tempo, que afeta principalmente o valor das premiações mais antigas. Com essa percepção, pode-se constatar que tal correção se revela necessária para uma modelagem adequada dos dados.

Com o aparato probabilístico apresentado e amplo uso da função “choose(n,r)” para calcular combinações, podemos construir duas tabelas que se encontram no site da CAIXA, a fim de subsidiar a análise. 

# Cálculo das chances e quantidades de prêmios a receber

s=c(6:15)
C0=4.5 # Custo base das apostas: valor da aposta simples de 6 números (mais barata)  
C=c()
d<-matrix(0,10,3)
q<-matrix(0,10,6)
for(i in 1:10){
  C[i]=C0*choose(i+5,6) # Vetor dos custos das apostas, formado por múltiplos do custo base
  for(j in 1:3){
    d[i,j]=round(choose(60,6)/(choose(i+5,7-j)*choose(55-i,j-1)),0) # Inverte probabilidades de acerto
    for(k in j:3){
      q[i,min(j+k,j*k)]=choose(7-j,7-k)*choose(i+j-2,k-1) # Prêmios a receber por faixa de acerto 
    }
  }
} # Contrução dos vetores que compõem as tabelas Probabilidade e Quantidade de Prêmios a Receber

chances <- data.frame(cbind(s,C,d))
colnames(chances)=c("Numeros_jogados","Valor_aposta","Probinv_sena","Probinv_quina","Probinv_quadra")

premios <- data.frame(cbind(s,q))
colnames(premios)=c("Numeros_jogados","Senas_faixa1","Quinas_faixa1","Quadras_faixa1",
"Quinas_faixa2","Quadras_faixa2","Quadras_faixa3")

Estas são as tabelas de probabilidades e prêmios no site da CAIXA:

Tabela de probabilidades Tabela de prêmios

Na sequência, nossa tabela construída, que combina as duas informações em uma, através de uma junção no R: 

# Matriz de probabilidades das apostas

probabilidades <- inner_join(chances,premios,"Numeros_jogados")
probabilidades <- probabilidades %>% mutate_at(vars(1,3:11),as.integer)
glimpse(probabilidades)
## Rows: 10
## Columns: 11
## $ Numeros_jogados <int> 6, 7, 8, 9, 10, 11, 12, 13, 14, 15
## $ Valor_aposta    <dbl> 4.5, 31.5, 126.0, 378.0, 945.0, 2079.0, 4158.0, 7722.0~
## $ Probinv_sena    <int> 50063860, 7151980, 1787995, 595998, 238399, 108363, 54~
## $ Probinv_quina   <int> 154518, 44981, 17192, 7791, 3973, 2211, 1317, 828, 544~
## $ Probinv_quadra  <int> 2332, 1038, 539, 312, 195, 129, 90, 65, 48, 37
## $ Senas_faixa1    <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1
## $ Quinas_faixa1   <int> 0, 6, 12, 18, 24, 30, 36, 42, 48, 54
## $ Quadras_faixa1  <int> 0, 0, 15, 45, 90, 150, 225, 315, 420, 540
## $ Quinas_faixa2   <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10
## $ Quadras_faixa2  <int> 0, 5, 15, 30, 50, 75, 105, 140, 180, 225
## $ Quadras_faixa3  <int> 1, 3, 6, 10, 15, 21, 28, 36, 45, 55

Modelagem do problema

Quando se faz um investimento, de qualquer natureza, é normal esperar um retorno, mas nunca sem risco. No caso da loteria, o valor investido é o desembolso com os cartões, enquanto o retorno é o prêmio desejado. Para o jogador, o custo da aposta é certo, fixado por quem promove os jogos de loteria, ao passo que o retorno é incerto, por ter uma probabilidade associada, a qual mensura o risco de ganho ou perda.

Essas noções podem ser resumidas no conceito chamado esperança estatística. Também denominada valor esperado, pode ser entendida como uma soma dos valores que uma determinada variável aleatória W de interesse pode assumir, ponderada pelas probabilidades de ocorrência de cada um desses valores. Essa descrição pode ser resumida na fórmula abaixo:

\[ E(W) = \sum_{j}w_jp(w_j) \]

De um modo geral, esse conjunto de índices j pode ser finito ou infinito. No contexto da Mega-Sena, a variável W em foco pode muito bem ser o prêmio a receber, enquanto os índices representam as diferentes faixas de premiação 1, 2 e 3, correspondendo respectivamente a Sena, Quina e Quadra. Além disso, como uma pessoa pode, em cada sorteio, jogar x cartões distribuídos alternativamente de diversas maneiras entre cada uma das 10 modalidades de aposta (6 a 15 números), é possível traduzir sua esperança de ganhar por meio da seguinte expressão matemática:

\[ E(W(x)) = E(W(x_1,\ldots,x_{10})) = \sum_{i=1}^{10}x_iE(W_i) \]

No último membro da expressão acima, cada nova parcela de esperança representa o prêmio esperado para cada modalidade de aposta. Como estamos tratando simultaneamente as possibilidades de acertar Sena, Quina ou Quadra, além de considerar o prêmio nulo para qualquer outra quantidade de acertos, os valores podem ser determinados da mesma forma que as matrizes de chances e ganhos: \[ E(W_i) = \sum_{j=1}^3w_{ij}p_{ij}, \ \ \ p_{ij}=\frac{C_{7-j}^{i+5}\times C_{j-1}^{55-i}}{C_{6}^{60}}, \ \ \ w_{ij}=\sum_{k=j}^3\big(C_{7-k}^{7-j}\times C_{k-1}^{i+j-2}\big)v_{k} \]

Apesar de não terem uma aparência muito agradável, todos os valores descritos acima estão tabelados, com exceção do vetor v que representa os prêmios associados a cada uma das 3 faixas de premiação. É aqui que entra o índice de inflação acumulada INPC, para corrigir todos os prêmios pagos ao valor atual e possibilitar o cálculo de sua média em cada faixa. Ao fazermos isso, obtemos o seguinte vetor de ganhos (equivalente a v):

# Cálculo do vetor de premiações médias referentes a cada faixa de aposta 

geral <- inner_join(resultados,inpc,"Periodo_referencial")

ganhos <- c(sum(geral$Rateio_faixa1*geral$INPC_acum_atual)/length(geral$Concurso),
     sum(geral$Rateio_faixa2*geral$INPC_acum_atual)/length(geral$Concurso),
     sum(geral$Rateio_faixa3*geral$INPC_acum_atual)/length(geral$Concurso))
print(ganhos)
## [1] 8661505.2058   47434.9483     770.8543

Será que ao olhar para os prêmios médios corrigidos no vetor de ganhos, ainda parece tão motivador jogar na Mega-Sena? Perceba, por exemplo, que a uma aposta simples com chance de 1 em 50.063.860 corresponde um prêmio médio aproximado pouco superior a 8 milhões de reais. Os outros dois valores também podem ser contrastados com as respectivas probabilidades na tabela correspondente. No entanto, como ainda é possível jogar mais cartões x e mais modalidades de aposta (de 1 a 10) em cada sorteio, começaremos a ver na sequência os efeitos dessa decisão.

Como a esperança estatística é linear, podemos construir a função esperança de ganhar na Mega-Sena no software R de duas maneiras: uma para visualização (E) e outra que atenda aos requisitos dos algoritmos de otimização (En). 

# Definição das funções de esperança

E <- function(x, P, g){
  mP <- as.matrix(P)
  vg <- as.vector(g)
  e <- c()
  p <- list(NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL)
  w <- list(NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL)
  for(i in 1:nrow(mP)){
    for(j in 1:length(vg)){
      p[[i]][j]=1/mP[i,j+2]
      w[[i]][j]=0
      for(k in j:length(vg)){
        w[[i]][j]=w[[i]][j]+mP[i,min(j+k+5,j*k+5)]*vg[k] 
      }
    }
    e[i]=sum(w[[i]]*p[[i]])
  }
  return(sum(e*x))
} # Valor esperado da premiação na Mega-Sena em termos do vetor x de cartões, da matriz de probabilidades e do vetor de ganhos

En <- function(x){
  e<-c()
  p<-list(NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL)
  w<-list(NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL)
  g<-c(8164791.288,44745.693,723.759)
  for(i in 1:10){
    for(j in 1:3){
      p[[i]][j]=choose(i+5,7-j)*choose(55-i,j-1)/choose(60,6)
      w[[i]][j]=0
      for(k in j:3){
        w[[i]][j]=w[[i]][j]+choose(7-j,7-k)*choose(i+j-2,k-1)*g[k] 
      }
    }
    e[i]=sum(w[[i]]*p[[i]])
  }
  return(-sum(e*x))
} # Função esperança negativa em termos do vetor x de cartões apostados para adequação aos algoritmos de minimização


# Lista de vetores canônicos para calcular esperança de prêmio em cada um dos 10 tipos de apostas (6 a 15 números)
b<-list()
u<-c()
for(i in 1:10){
  for(j in 1:10){
    if(j==i){
      u[j]<-1
    }else{
      u[j]<-0
    }
  }
  b[[i]]<-u
}

É possível vislumbrar o comportamento dessa função em cada modalidade de aposta, avaliando-a na base canônica do espaço euclidiano de 10 dimensões, construída numa lista de vetores que denominamos b. Obviamente não dá para visualizar o gráfico completo de E, o que apenas ocorreria em um espaço com 11 dimensões. Isso lembra uma certa teoria… 

# Visualizações da função E avaliada para cada tipo de aposta

for(i in 1:10){
  plot(i+5,E(b[[i]],probabilidades,ganhos)/1000,main="Valor esperado em cada tipo de aposta",xlab="Números jogados",xlim=c(6,15),ylab="Esperança do prêmio (mil R$)",ylim=c(0,4),pch=19,col="blue")
  par(new = T)
}

Ao contemplar o gráfico da esperança, parece muito promissor jogar mais números, não? Entretanto, nada vem de graça. Em um problema de otimização, essa máxima equivale à presença de restrições. No nosso contexto específico, faz sentido pensar em pelo menos três limitações básicas.

Primeiramente, devemos assumir que o jogador possui uma quantia máxima R para gastar com as apostas, cujo preço em cada modalidade está expresso no vetor C de custos das apostas, determinado anteriormente. Ademais, também é razoável estabelecer uma quantia máxima N de cartões a jogar em um sorteio. Por fim, ainda há limites para as variáveis em decorrência da sua natureza discreta. Tais restrições podem ser expressas matematicamente como segue:

Na sequência, mostramos um esboço de tal conjunto de restrições em 3 dos 45 planos de projeção possíveis, para N = 20.000 e R = R$ 500.000,00. O conjunto de restrições, de acordo com a ordem introduzida acima, possui as cores laranja, lilás e preta.

Restrições da Mega-Sena

Visualmente dá para notar, se ainda não ficou claro, que todas as restrições são lineares, bem como a própria função objetivo. Além disso, embora o conjunto viável tenha sido representado de forma contínua nos gráficos, ele é formado por pontos isolados (conjunto discreto), uma vez que o problema apenas tem sentido no domínio dos números inteiros não negativos.

Dessa forma, fica determinado um problema de otimização linear inteira (PLI) para expressar a busca pela melhor estratégia do jogador na Mega-Sena:

\[ \mbox{ maximizar } \ E(x) = \sum_{i=1}^{10}{x_i}E(W_i) \] \[ \mbox{ sujeito a } \ \sum_{i=1}^{10}x_iC_i - R \le 0, \] \[ \ \ \ \ \ \sum_{i=1}^{10}x_i - N \le 0, \] \[ \ \ \ \ \ 0 \le x_i \le \frac{C_{6}^{60}}{C_{6}^{i+5}},\ x_i \in \mathbb{Z}. \]

Pode-se verificar que o conjunto de restrições delimitado pelas funções lineares é compacto e a função objetivo E é claramente contínua, donde atendem ao Teorema de Weierstrass. Este garante a existência de valor máximo e valor mínimo para E (e, por conseguinte, para a função En).

Abstraindo o que seja um conjunto compacto ou uma função contínua, o importante é saber que, matematicamente, o problema tem solução (não necessariamente única). Outro detalhe a enfatizar é que, a rigor, a otimização ocorre em um conjunto discreto com um número finito de elementos e não em uma região contínua do espaço euclidiano. Esse fato não altera as propriedades fundamentais que garantem a existência de solução; entretanto, por mais contraintuitivo que pareça, a restrição ao domínio dos inteiros torna o problema mais difícil (principalmente para valores elevados das constantes R e N).

Conforme a natureza peculiar do problema, o método utilizado será a Programação Linear Inteira, através da função “lp” do pacote “lpSolve” do R.

A fim de fazer testes e comparações, verificaremos três cenários:

  1. Limites orçamentário e de cartões baixos: R1 = R$ 300,00 e N1=20.
  2. Limites orçamentário e de cartões médios: R2 = R$ 1.000,00 e N2=100.
  3. Limites orçamentário e de cartões altos: R3 = R$ 25.000,00 e N3=500.

Na sequência, mostramos as constantes globais, aplicadas a todos os cenários. 

# Limites globais

C <- probabilidades$Valor_aposta # Vetor de custos   
apmax <- probabilidades$Probinv_sena # Limites combinatórios  
A <- rbind(C,rep(1,10),diag(1,10,10)) # Matriz de restrições


Tendências e loterias - uma imagem diz mais que mil contas!

Cenário estratégico 1

Mote: Menos é mais…

Faço minha fezinha mesmo!

# Otimização das apostas no cenário 1

R1=300
N1= 20
iL1<-Sys.time()
S1<-lp(direction = "min", objective.in = c(En(b[[1]]),En(b[[2]]),En(b[[3]]),En(b[[4]]),En(b[[5]]),
En(b[[6]]),En(b[[7]]),En(b[[8]]),En(b[[9]]),En(b[[10]])), const.mat = A, const.dir =
c("<=","<=","<=","<=","<=","<=","<=","<=","<=","<=","<=","<="), const.rhs = c(R1, N1, apmax), all.int = T)
fL1<-Sys.time()
y1<-S1$solution 
print(y1) # Solução encontrada: (10,8,0,0,0,0,0,0,0,0) -> 18 cartões

10, 8, 0, 0, 0, 0, 0, 0, 0, 0


Tabela de avaliação dos resultados para o primeiro cenário: 

# Avaliação dos retornos no cenário 1: 

h1 <- fL1-iL1 # Tempo médio de execução do algoritmo
M1 <- -En(y1) # Valor da esperança de ganho otimizada

# Avaliação dos custos no cenário 1: 

C1 <- sum(C*y1) # Valor do custo final
I1 <- -En(y1)/sum(C*y1) # Recuperação do investimento
Métricas Monetárias Cartões Execução (média de 10x)
A) Retornos ou Saídas R$ 50,36 18 -
B) Custos ou Limites R$ 297,00 20 0,006402683 seg.
C) Relativo: A/B 16,9552% 90% -

O ponto ótimo para o qual o algoritmo converge neste cenário aproxima-se tanto da fronteira de gasto (99% da restrição orçamentária) quanto da fronteira do número de cartões. Observe que o algoritmo tem relativamente poucas condições a testar, uma vez que o limite de orçamento do jogador não lhe permite fazer cartões com mais de 8 números. Assim, o mote parece estar de acordo com a estratégia possível do apostador, embora o algoritmo entregue um retorno esperado ótimo relativo ao custo não muito animador.

Cenário estratégico 2

Mote: O universo tende ao equilíbrio.

Pensando na confluência cósmica da numerologia.

# Otimização das apostas no cenário 2

R2=1000
N2=100
iL2<-Sys.time()
S2<-lp(direction = "min", objective.in = c(En(b[[1]]),En(b[[2]]),En(b[[3]]),En(b[[4]]),En(b[[5]]),
En(b[[6]]),En(b[[7]]),En(b[[8]]),En(b[[9]]),En(b[[10]])), const.mat = A, const.dir =
c("<=","<=","<=","<=","<=","<=","<=","<=","<=","<=","<=","<="),const.rhs = c(R2, N2, apmax), all.int = T)
fL2<-Sys.time()
y2<-S2$solution
print(y2) # Solução encontrada: (26,8,2,1,0,0,0,0,0,0) -> 37 cartões

26, 8, 2, 1, 0, 0, 0, 0, 0, 0


Tabela de avaliação dos resultados para o segundo cenário: 

# Avaliação dos retornos no cenário 2: 

h2 <- fL2-iL2 # Tempo médio de execução do algoritmo
M2 <- -En(y2) # Valor da esperança de ganho otimizada: R$ 169,38

# Avaliação dos custos no cenário 2:

C2 <- sum(C*y2) # Valor do custo final
I2 <- -En(y2)/sum(C*y2) # Recuperação do investimento
Métricas Monetárias Cartões Execução (média de 10x)
A) Retornos ou Saídas R$ 169,38 37 -
B) Custos ou Limites R$ 999,00 100 0,01278162 seg.
C) Relativo: A/B 16,9552% 37% -

A convergência do algoritmo neste cenário dá-se para um jogo 99,9% próximo da fronteira de gasto e com um número de cartões relativamente baixo em relação ao limite do jogador, além de tender ligeiramente para mais apostas com menos números. Porém, o que mais chama a atenção é o retorno ótimo relativo ao custo igual ao do cenário 1. Será que esse resultado sugere algum tipo de equilíbrio lotérico?

Cenário estratégico 3

Mote: Viva Las Vegas!

Vou arriscar tuudooo!!!

# Otimização das apostas no cenário 3

R3=25000
N3=500

iL3<-Sys.time()
S3<-lp(direction = "min", objective.in = c(En(b[[1]]),En(b[[2]]),En(b[[3]]),En(b[[4]]),En(b[[5]]),
En(b[[6]]),En(b[[7]]),En(b[[8]]),En(b[[9]]),En(b[[10]])), const.mat = A, const.dir =
c("<=","<=","<=","<=","<=","<=","<=","<=","<=","<=","<=","<="), const.rhs = c(R3, N3, apmax), all.int = T)
fL3<-Sys.time()
y3<-S3$solution 
print(y3) # Solução encontrada: (25,5,5,1,2,2,1,0,1,0) -> 42 cartões

25, 5, 5, 1, 2, 2, 1, 0, 1, 0

Tabela de avaliação dos resultados para o terceiro cenário: 

# Avaliação dos retornos no cenário 3: 

h3 <- fL3-iL3 # Tempo de execução do algoritmo
M3 <- -En(y3) # Valor da esperança otimizada: R$ 4.238,37

# Avaliação dos custos no cenário 3:

C3 <- sum(C*y3) # Valor do custo final: R$ 24.997,50
I3 <- -En(y3)/sum(C*y3)
Métricas Monetárias Cartões Execução (média de 10x)
A) Retornos/Saídas R$ 4.238,37 42 -
B) Custos/Limites R$ 24.997,50 500 317,71002 seg.
C) Relativo: A/B 16,9552% 8,4% -

Desta vez, o algoritmo convergiu para o ótimo, atingindo o ponto mais próximo da fronteira de gasto entre os três cenários (99,99% da restrição orçamentária) e a menor relação entre o número de cartões ótimo e o limite do apostador. Novamente, destaca-se no ponto ótimo a prevalência de cartões com menos números, mesmo num cenário em que o jogador pode apostar em todas as modalidades. É marcante também a grande elevação do tempo de execução do algoritmo, apontando para uma provável tendência exponencial de crescimento desse tempo com o aumento das constantes nas restrições. O mais interessante, porém, está no valor da esperança retorno ótima relativa ao custo: tão baixo quanto e idêntico ao obtido nos cenários anteriores. Coincidência?

Conclusões e extensões possíveis

Em linhas gerais, pode-se afirmar que esta análise atingiu seu objetivo, pois é possível responder às perguntas inicialmente formuladas, na forma de algumas conclusões relevantes, que são elencadas abaixo:

  1. Pode-se dizer que o coeficiente de recuperação do investimento associado às apostas na Mega-Sena, o qual relaciona o valor ótimo da esperança na aposta ao seu respectivo custo, é aproximadamente constante. Seu valor equivale a menos de 17% e revela um limite baixo para a esperança.
  2. Todos os resultados ótimos ficaram muito próximos da fronteira do gasto (99% ou mais), o que aponta para uma dependência direta muito forte entre o aumento do gasto com as apostas e o aumento da chance de ganhar, contribuindo para corroborar o item anterior.
  3. As soluções ótimas do algoritmo de PLI parecem distribuir os jogos de maneira a priorizar a realização de mais apostas com menos números, o que pode representar um mecanismo implícito para dispersar mais os números entre os cartões.


As conclusões 1 e 2 juntas sugerem que, em nenhum cenário, o jogador supera os limites intransponíveis da esperança relativa de 17%, a qual corresponderia à procurada relação retorno x custo das apostas. Isso responde à questão 2.

Ademais, baseado na conclusão 3, existe uma melhor estratégia do jogador no sentido de ajustar o nível de variabilidade das apostas com restrições de qualidade, tais como o controle do número de cartões a apostar, através de um algoritmo de PLI. Isso responde à questão 1 (pelo menos parcialmente).

Como extensões a investigar, poder-se-ia fazer uma conjectura mais geral sobre o coeficiente de recuperação do investimento na Mega-Sena (ou outras loterias): esse valor seria uma constante para qualquer padrão de jogo adotado na mesma loteria, seja ele ótimo ou não, inclusive em bolões?

Por fim, caberia ainda o uso de análise de Séries Temporais e técnicas de Machine Learning para avaliar o possível desenho de uma estratégia de jogo otimizada ao longo do tempo.