Modelo de Bornholdt para simulação de mercados financeiros artificiais

Grupo: Leonardo Barcelos, Luana Bianchi e Rubens Borrasca

Em 2002, Stefan Bornholdt apresentou um modelo de spins motivado pela dinâmica de agentes no mercado financeiro, que retornava observações observadas em situações reais, como distribuição de retorno de caudas pesadas e clusterização de volatilidade.^[1]

O objetivo deste trabalho é, então, apresentar o Modelo de Bornholdt como alternativa para a construção de um mercado financeiro artificial com características observadas em mercados reais. Dentre os principais resultados obtidos, destaca-se a formação de estados metaestáveis, distribuição de retornos de cauda pesada, e clusters de volatilidade.

Por fim, também é apresentado uma alternativa para o uso do Modelo de Bornholdt, analisando a opinião dos agentes da rede, e como estas influenciam no resultado final.

Introdução

Para estudar os fenômenos físicos que nos rodeiam, usamos modelos matemáticos para entender como sistemas evoluem com o tempo. No entanto, esses modelos muitas vezes podem ser base para estudos além do sistema físico de interesse. Para um sistema de spins, um dos modelos mais simples que leva em conta apenas a interação de cada spin com seus primeiros vizinhos, é o Modelo de Ising^[2]. O hamiltoniano que o descreve é defindo pela seguinte expressão:

 ${\mathcal {H}}=\sum _{\langle i,j\rangle }J_{ij}S_{i}^{z}S_{j}^{z},$

onde $J_{ij}$ é o termo que define a magnitude de interação entre os spins i e j e $S^{z}$ a componente Z do spin. Além disso $\langle i,j\rangle$ denota que estamos somando somente sobre os primeiros vizinhos. Para um valor de $J_{ij}>0$ , a interação ferromagnética é favorecida, e os spins vizinhos tendem a se alinhar conjuntamente, formando domínios magnéticos. No entanto, para $J_{ij}<0$ a interação antiferromagnética é favorecida, e os spins tendem a se "anti-alinhar".

Simulação de Mercados de Ações e o Modelo de Bornholdt

Podemos fazer um paralelo interessante entre os spins e a sua interação com seus vizinhos com traders em mercados finaceiros e a suas interações com outros traders. No caso do Modelo de Ising, podemos interpretar cada spin como um trader, ou agente, em um mercado de ações e a direção da componente Z desse spin como a sua estratégia de compra ou venda. Então neste caso mais simples, para $J_{ij}>0$ , os agentes tenderão a adotar a mesma estratégia de mercado, e isso não descreve bem a realidade dos mercados. Para o caso em que $J_{ij}<0$ , os agentes irão adotar estratégias contrárias aos seus vizinhos, e isso também não nos levará a uma representação realista do mercado financeiro. Para que o modelo possa se aproximar da realidade, devemos incrementar o hamiltoniano do sistema.

Um modelo bastante próximo da realidade dos mercados de ações é o proposto por Lux e Marchesi ^[3] ^[4] que classificam os agentes em duas estratégias: os fundamentalistas e os chartistas. Este modelo reproduz muitas das propriedades observadas em mercados reais, como lei de distribuição de potências (cauda pesada) dos retornos dos preços e uma alta correlação da volatilidade de preços. Além disso, no nível de estratégias, apresenta o fenômeno em que o número de chartistas se correlacionam com fases de grande volatilidade, como visto em mercados reais. Este modelo tem grande sucesso, no entanto sua complexidade é alta. Stefan Bornholdt então propôs um modelo de máxima simplicidade, baseado no modelo de Ising, desenhado para simular a dinâmica de expectativas em sistemas de muitos agentes.

No modelo de Bornholdt, há pelo menos duas forças conflitantes vistas em ações econômicas:

O quê os vizinhos fazem: comumente associado à ação de chartistas;
O quê a minoria faz: comportamento associado a fundamentalistas, traders com maior conhecimento sobre valores fundamentais dos mercados de ações.

Mais detalhes sobre como a opinião de um agente influencia no sistema estão na seção de opinião dos agentes

Neste modelo, essas duas interações conflitantes são combinadas: as interações entre vizinhos são representadas pelo modelo de Ising mais simples; um acoplamento à minoria como um observável global é introduzido por um acoplamento à magnetização global do sistema de spins. Assim, o hamiltoniano incrementado fica com a seguinte expressão:

 ${\mathcal {H}}=\sum _{\langle i,j\rangle }J_{ij}S_{i}S_{j}-\sum _{i}\alpha C_{i}{\frac {1}{N}}\sum _{j}S_{j},$

onde $\alpha$ é o termo de acomplamento de cada spin com a magnetização da rede $(\alpha >0)$ e $C_{i}$ representa a estratégia com relação à magnetização da rede.

Por simplicidade, assume-se que cada spin é atualizado com uma dinâmica de banho térmico de acordo com:

 $S_{i}(t+1)={\begin{cases}+1,&{\text{com }}p=1/[1+exp(-2\beta h_{i}(t))]\\-1,&{\text{com }}1-p\end{cases}},$

onde $\beta ={\frac {1}{k_{B}T}}$ . Neste trabalho, por simplicidade assumiu-se que $k_{B}=1$ .

Considerando um modelo com $i=1,...,N$ spins, com as orientações $S_{i}(t)=\pm 1$ , a dinâmica dos spins dependerá do campo local $h_{i}(t)$ :

 $h_{i}(t)=\sum _{j}^{N}J_{ij}S_{j}-\alpha C_{i}(t){\frac {1}{N}}\sum _{j}^{N}S_{j}(t).$

Neste trabalho, o termo de interação é tomado como uma constante $J_{ij}=1$ para vizinhos e $J_{ij}=0$ para os demais spins, e ele é responsável pela indução de ordem ferromagnética local.

Cenários de Simulação

Vários cenários podem ser analisados, dependendo do valor e da dinâmica das estratégias. Consideremos o caso mais simples, em que $C_{i}(t)=1\forall i,t$ . Cada agente, apesar do acoplamento ferromagnético local com os primeiros vizinhos, tem um acoplamento antiferromagnético com a magnetização. Esta dinâmica corresponde a traders que em adição a um nível básico de ferromagnetismo, ou seja, tendem a adotar estratégias parecidas com seus vizinhos, também têm o desejo de se juntar a minoria global, por exemplo a fim de investir em possíveis ganhos futuros. Portanto, traders com $C_{i}(t)=1$ podem ser chamados de fundamentalistas. Se todos os agentes adotam esta estratégia, a dinâmica global do sistema tenderá rápidamente a um estado de magnetização quase zero, mesmo para temperaturas abaixo da temperatura crítica $T_{c}$ . Podemos interpretar a temperatura nesse sistema como se fosse a velocidade com que as transações acontecem: para uma temperatura baixa, os agentes trocariam suas estratégias mais devagar, por exemplo. A temperatura no caso deste modelo, vai ditar a distribuição de estratégias: para uma temperatura baixa, as estratégias associadas ao spin up dominarão; já para uma temperatura mais elevada, as estratégias associadas ao spin down predominarão.

Um cenário mais interessante, é permitir que os agentes possam adotar duas estratégias diferentes, sendo possível agora adotar $C_{i}(t)=-1$ , que corresponde a um acoplamento ferromagnético com a magnetização global. Essa estratégia é chamada de chartista, visto que os agentes tendem a seguir a opinião da maioria dos traders.

Por fim, podemos definir regras para a transição entre as duas estratégias estratégias, onde cada trader tenderá a adotar uma estratégia ótima. Consideremos o cenário mais simples para trocas de estratégia: um agente no grupo majoritário frequentemente tenderá a mudar para o grupo minoritário, por exemplo para apostar em um comodity que ainda não está na moda (e possivelmente escapar de um crash do seu bem mais popular no momento). Por outro lado, um agente que se encontra no grupo minoritário (portanto esperando retornos futuros) pode não estar satisfeito com seus retornos atuais. Em resumo agentes nos grupos majoritários sempre irão escolher a estratégia $C_{i}(t)=1$ , enquanto a minoria irá escolher $C_{i}(t)=-1$ . Cada agente escolhe uma estratégia arriscada a fim de aumentar seus retornos. A dinâmica da transição de estratégias é dada pela expressão:

 $C_{i}(t+1)=-C_{i}(t),{\text{ se }}\alpha C_{i}(t)\sum _{j=1}^{N}S_{j}(t)<0.$

Em todos os casos, a simulação foi feita usando a dinâmica de banho térmico já apresentada, juntamente com Monte Carlo. Para cada passo (sweep) de MC, sorteia-se várias posições no sistema e aplica-se a dinâmica de heat-bath sobre cada uma, calculando, após isso, a orientação do spin.

Alguns conceitos importantes

Retornos ^[5]

Quando se trata de sistemas financeiros, os estudos se concentram mais no retorno dos ativos do que no preço em si, pois a série temporal dos retornos tem propriedades estatísticas mais interessantes que a série dos preços.

Sendo P(t) o preço de um ativo financeiro no instante t, e P(t-1) o preço do ativo no instante (t-1), o retorno linear do ativo é:

 $r(t)={\frac {P(t)-P(t-1)}{P(t-1)}}$

Reescrevendo esta equação, obtemos que:

 $r(t)+1={\frac {P(t)}{P(t-1)}}$

Aplicando a função logarítmica em ambos os lados da equação, e considerando que:

 $ln(x+1)\approx x$

obtêm-se o retorno logarítmico, que é mais indicado quando se têm ativos voláteis, que possuem uma variação muito alta:

 $r(t)=ln\left({\frac {P(t)}{P(t-1)}}\right)$

Considerando que neste estudo serão comparados retornos de diferentes índices, e também os retornos obtidos através das simulações com o modelo de Bornholdt, é importante normalizar os retornos, para que fiquem na mesma escala:

 ${\bar {r}}_{n}={\frac {r_{n}-\langle r\rangle }{\sigma _{r}}}$

em que ${\textstyle \sigma _{r}}$ é o desvio padrão da serie de retornos e ${\textstyle \langle r\rangle }$ a média.

Distribuição dos Retornos ^[6]

Quando se tem um volume considerável de dados é possível obter a distribuição probabilística deles. Para isso pode-se utilizar a estimação de densidade de Kernel (KDE)^[6]. Ao observar uma pequena janela de tamanho 2h em torno de um ponto em análise, pode-se dizer que:

 $P(x-h\leq x\leq x+h)\approx {\frac {1}{N}}\sum _{n}^{N}W(u)$

sendo $W(u)$ uma função kernel e $u$ uma variável tal que:

 $u={\frac {x-x_{n}}{2h}}$

Para este estudo utilizou-se um kernel gaussiano:

 $W(u)={\frac {1}{\sqrt {2\pi }}}e^{-{\frac {1}{2}}u^{2}}$

Este método foi aplicado para as séries de retorno para obter a distribuição deles, utilizando $h=0.5\sigma$ .

Volatilidade

Exemplo de uso da volatilidade como série temporal de um retorno arbitrário.

Uma forma de calcular a volatilidade da série temporal de retornos ao longo do tempo é elevar ao quadrado os valores da série. Deste modo pode-se obter uma variável como a que está ilustrada na figura ao lado.

O interessante em estudar volatilidade de retornos financeiros é que essa variável reflete o quão imprevisível é um determinado ativo. Uma ação com alta volatilidade tende a ter um risco maior de investimento, ao passo que ações com baixa volatilidade geralmente retornam riscos menores, pois seu comportamento acaba sendo mais previsível.

Um fato estilizado financeiro é que a volatilidade das séries temporais de retorno apresentam comportamento sazonal por natureza. ^[1] ^[7] ^[4] Há períodos de alta volatilidade, seguidos por períodos com baixa volatilidade, que então são novamente seguidos por alta volatilidade, e assim adiante. E uma forma de mensurar isto é verificando a presença de clusters de autocorrelação na volatilidade de retornos. Isto é, através da análise da autocorrelação da volatilidade, encontrar bolhas que indiquem as fases destes comportamentos.

Para obter a auto correlação o Teorema de Wiener-Khinchin ^[8] ^[9] foi utilizado, de forma que:

 $A(\tau )={\mathcal {F}}\{|S(\omega )|^{2}\}$

onde $S(\omega )$ é a Transformada de Fourier do quadrado dos retornos.

Simulações

Variação do tamanho da grade

Um dos objetivos deste estudo é verificar qual o tamanho de grade gera um resultado que melhor simula um índice financeiro. A escolha de índices financeiros para comparação, ao invés de ações ou commodities, é por causa da instabilidade que ativos financeiros separados têm. Imaginando um cenário fictício onde uma empresa A vende sorvetes enquanto a empresa B vende chocolate quente. A tendência é que haja uma sazonalidade nos 2 ativos, de tal forma que, no verão, as ações de A subam enquanto B desce, e vice-versa no inverno. Porém, enquanto os 2 ativos estão se movimentando, o movimento do mercado representado por estas 2 ações permanece aproximadamente constante. Como índices financeiros são, na verdade, médias de um conjunto grande e diverso de ações presentes no mercado, são melhores para a comparação com o modelo.

Foram escolhidos 4 diferentes tamanhos de grade:

16x16: 256 agentes
32x32: 1024 agentes
50x50: 2500 agentes
100x100: 10000 agentes

Cada uma das simulações se deu sob as seguintes condições:

$\alpha =4$
$J=1$
$T=1,5$
Todos os agentes são racionais, ou seja, eles podem trocar de opinião ao longo da simulação, assumindo $C_{i}=-1$ em alguns momentos e $C_{i}=1$ em outros, de acordo com o grupo que o agente está inserido: maioria ou minoria.

Abaixo é possível observar fotografias do sistema em diferentes instantes para cada uma das simulações:

Situação do sistema de 256 agentes (grade 16 x 16) em 9 instantes diferentes.	Situação do sistema de 1024 agentes (grade 32 x 32) em 9 instantes diferentes.	Situação do sistema de 2500 agentes (grade 50 x 50) em 9 instantes diferentes.	Situação do sistema de 10000 agentes (grade 100 x 100) em 9 instantes diferentes.

Como saída da simulação é obtida uma série temporal dos valores de magnetização, que neste modelo significam o preço de um ativo presente no mercado. Utilizando a série de preços (magnetização), foi obtido, então, a série dos retornos logarítmicos normalizados, através do método explicado na seção dos Retornos. Para comparar o modelo com dados reais, obteve-se as series temporais do preço de índices do mercado financeiro: Ibovespa, S&P 500, Dow Jones e NASDAQ. Esses dados foram obtidos com a biblioteca Pandas Datareader do Python, que permite obter dados de páginas da internet como a Yahoo Finance, que possui a série temporal dos preços de vários ativos financeiros. Assim como nas simulações, a partir da série de preços dos índices foram calculadas a série temporal dos retornos logarítmicos normalizados.

Na figura abaixo pode-se observar os retornos de cada simulação e índice financeiro:

Retornos logarítmicos normalizados obtidos para as simulações com diferentes números de agentes (tamanho da grade) na primeira coluna, e na segunda coluna os retornos logarítmicos normalizados para os índices financeiros Ibovespa, S&P 500, Dow Jones e NASDAQ

Alguns pontos interessantes a se observar:

os retornos dos índices S&P 500, Dow Jones e NASDAQ são muito parecidos, isso de deve ao fato de que todos são índices de mercados dos Estados Unidos, enquanto o Ibovespa é um índice de mercado brasileiro;
os retornos das simulações vão variando mais tempo perto do zero conforme o número de agentes cresce;
conforme o número de agentes no sistema aumenta, diminui-se o número de transições agudas de retorno, representadas pelos picos no gráfico. Isso indica um mercado mais estável.
o objeto da simulação não é retornar uma cópia do observado nos dados reais, mas sim um cenário parecido. O intuito é criar um novo mercado, artificial, com características parecidas ao observado nos mercados reais.

Para ter uma comparação melhor a fim de entender que tamanho de grade simula melhor um mercado financeiro, é importante ver a distribuição dos retornos e a auto correlação das volatilidades. Com base no conteúdo da seção Distribuição dos Retornos foi gerado uma curva de densidade de probabilidade para cada série de retorno, que podem ser observadas na figura abaixo:

Distribuição dos retornos das simulações (linha cheia) e dos índices financeiros (linha tracejada). Percebe-se que a distribuição dos retornos das simulações que mais se aproximam da distribuição dos índices financeiros são as que possuem um menor número de agentes.

Desta figura pode-se também observar alguns pontos:

quanto maior o número de agentes, mais longe a distribuição dos retornos das simulações ficam da distribuição dos retornos dos índices de mercado;
as caudas das distribuições dos retornos das simulações vão ficando mais pesadas conforme o número de agentes diminui, o que indica mais uma vez que as simulações com menor número de agentes possuem distribuição dos retornos mais parecidas com a do mercado financeiro, que também possui uma cauda pesada;
para números muito grandes de agentes, há uma distribuição de retorno de cauda mais leve, o que implica em muito menos situações de alto retorno, um mercado mais estável e consequentemente irreal.

Na figura abaixo estão presentes as auto correlações dos retornos quadrados, ou seja, da volatilidade dos retornos:

Auto correlação das volatilidades dos retornos das simulações (primeira coluna) e do retorno dos índices financeiros (segunda coluna). Em ambas as colunas, a escala dos eixos y e x dos gráficos é logarítmica.

Desta imagem são notados alguns aspectos como:

o comportamento da auto correlação das volatilidades das simulações com mais agentes são mais parecidas com as dos índices financeiros;
embora a simulação com 256 agentes tenha a distribuição dos retornos que mais se aproximou, analisando a sua auto correlação da volatilidade se observa que para este número de agentes os retornos são muito voláteis, e portanto não descrevem tão bem um sistema financeiro;
mesmo que as simulações com maior número de agentes apresente a auto correlação das volatilidades mais similar as do mercado financeiro, a partir da simulação com 1024 agentes, nota-se a presença de clusters de volatilidade, comprovando o fato estilizado de que existe memória do retorno quadrado.

Considerando os pontos apresentados, percebe-se que grades de um tamanho grande não descrevem o mercado financeiro da melhor forma, e da mesma forma, mesmo que aparente descrever bem devido a distribuição de retornos semelhantes a do mercado financeiro, a simulação de 256 agentes é muito volátil se comparada com um mercado real. Com isto, considera-se que a grade de tamanho 32x32 melhor representa um mercado de ações, pois possui uma distribuição dos retornos que não difere tanto da distribuição dos retornos dos índices financeiros, e ainda apresenta o fato estilizado da memória da volatilidade.

Variação da opinião dos agentes

O parâmetro $C_{i}$ , como dito anteriormente, indica a opinião de um agente $i$ presente na rede. Esta opinião está relacionada a seguir (ou não) o comportamento da maioria dos outros agentes presentes no sistema.

Por consequência, pode-se resumir as possíveis opiniões que um agente tem na rede em 3 ramos diferentes:

Opinião 1: $C_{i}=+1$ . Agente que opõe sua opinião à da maioria dos outros agentes no sistema, chamado de fundamentalista. Recebem este nome pois se apoiam no princípio fundamental de oferta e demanda da economia, o qual diz que, com demanda maior, a oferta é menor, enquanto que para demandas menores, a oferta é maior. Logo, opinando diferente da maioria, garante-se um maior retorno.
Opinião 2: $C_{i}=-1$ . Agente que iguala sua opinião à da maioria dos outros agentes no sistema, chamado de chartista. Recebe este nome do termo em inglês chart (gráfico). São agentes que sempre analisam, através de gráficos, as ações em alta, para comprá-las, e as em baixa, para vendê-las. Por isso, seguem a maioria do sistema.
Opinião 3: $C_{i}=c\in \{-1;+1\}$ . Agente que não possui estratégia, e a cada passo da simulação joga aleatoriamente do lado dos chartistas ou fundamentalistas, sem raciocínio. É chamado de completamente irracional.

Abaixo, estão gráficos que mostram o comportamento do sistema composto inteiramente por cada tipo de opinião considerada na rede.

Simulação com 100% da rede contendo agentes de opinião 1	Simulação com 100% da rede contendo agentes de opinião 2	Simulação com 100% da rede contendo agentes de opinião 3

Para o sistema composto apenas pela opinião 1, nota-se uma volatilidade muito grande na série temporal. Como os agentes de opinião 1 sempre jogam contra a maioria, quando o sistema atinge um estado de preços grande o suficiente, por exemplo, ocorre uma "debandada" dos agentes, que passam a opinar contra. Desta forma, o preço abaixa rapidamente, até chegar num estado suficientemente baixo, onde os agentes passam a opinar a favor. Não é um sistema muito realista pois os retornos não oscilam tão rapidamente assim.

Já no sistema composto apenas da opinião 2, há uma estabilidade quase que instantânea. Como todos os agentes seguem a maioria, uma vez atingido um determinado número mínimo de agentes com a mesma opinião, todos os outros passam a segui-la. E como não há agentes que discordem dessa opinião na rede, o sistema se mantem neste estado ad eternum. É um sistema bem longe do realista, pois há praticamente 0 risco.

Por último, o sistema composto da opinião 3 é, dentre os 3, o que mais se parece com a realidade. Isso porque, mesmo que de forma desorganizada e sem sentido lógico, os agentes acabam tendo opiniões diferentes, o que resulta em ganhadores e perdedores.

A imagem abaixo mostra as distribuições de retorno para cada um dos 3 casos mencionados:

Distribuição de retornos para 100% da rede contendo agentes de opinião 1	Distribuição de retornos para 100% da rede contendo agentes de opinião 2	Distribuição de retornos para 100% da rede contendo agentes de opinião 3

Enquanto que, para a rede composta apenas de agentes com opinião 1, temos caudas tão pesadas que apresentam sub-picos, na rede composta de opinião 2 a cauda é leve, e basicamente todos os retornos estão situados bem próximos à média da distribuição. Por último, para o caso 3, a cauda é pesada, porém ainda há concentração grande de retornos muito próximos ao centro da distribuição.

Com o intuito de construir uma rede de opiniões mistas (como é em casos de mercados reais), foram testadas várias combinações de porcentagens de opiniões na rede, até se encontrar uma que reproduzisse resultados muito semelhantes à de dados reais.

Utilizando aproximadamente 80% de agentes com opinião 3, 10% de opinião 1 e 5% de opinião 2, obteve-se o resultado ilustrado nos gráficos abaixo:

Simulação com 15% da rede contendo agentes de opinião 1, 5% de opinião 2 e 80% de opinião 3	Distribuição de retornos com 15% da rede contendo agentes de opinião 1, 5% de opinião 2 e 80% de opinião 3

Alguns pontos importantes:

provavelmente, se fossem feitos mais testes de combinações de probabilidade, chegaria-se a um resultado que encaixaria muito bem com alguma das curvas de dados reais. Porém, o intuito do trabalho, como dito anteriormente, não é replicar exatamente o que acontece com o mercado X ou Y, mas sim produzir um mercado artificial com características muito parecidas às observadas em mercados reais.
os retornos apresentam transições de estados de alta e baixa muito mais agudas do que casos reais. É quase como se não houvesse meio termo, ou o agente ganha muito, ou perde muito, ou não ganha nada. Uma limitação do modelo.
a distribuição de retornos apresenta cauda pesada, como observado em mercados reais.

A simulação acima indica que, para melhor simular mercados reais, a maioria dos agentes devem atuar de maneira irracional, enquanto que uma parcela de aproximadamente 10% atua de acordo com a minoria e aproximadamente 5% com a maioria.

Pode-se, também, permitir que os agentes troquem de opinião de forma minimamente racional ao longo da simulação. Deste modo, sempre que um agente "errar" no palpite (i.e. trocar a sua orientação de spin), ele escolhe outra opinião para seguir, até que erre novamente. A escolha da opinião é feita de forma aleatória pelo agente (portanto, não é 100% racional). O resultado para diferentes condições iniciais encontra-se nas figuras abaixo:

Porcentagem de composição da rede para cada opinião ao longo da simulação, com condição inicial igual para as 3 opiniões.	Porcentagem de composição da rede para cada opinião ao longo da simulação, com condição inicial de 100% opinião 3.	Porcentagem de composição da rede para cada opinião ao longo da simulação, com condição inicial de 100% opinião 1.	Porcentagem de composição da rede para cada opinião ao longo da simulação, com condição inicial de 100% opinião 2.

Observando a figura acima, nota-se a tendência do sistema em colapsar para um estado onde a maioria dos agentes têm opinião 2, até mesmo quando inicia-se com 0 agentes partilhando desta opinião.

Relembrando, a opinião 2 refere-se aos agentes que têm opinião ajustada pela maioria. Em todas as simulações feitas, após um certo tempo transcorrido da simulação, 100% da rede é composta por esta opinião, o que implica em um estado de consenso da rede quanto à opinião de compra ou venda sobre um ativo. Faz muito sentido que a tendência natural do sistema seja colapsar para esta opinião, porque nela todos os agentes saem sem perder dinheiro da negociação (em compensação, também ganham um retorno ínfimo).

O problema é que este resultado não pode ser observado na prática, pois o modelo de Bornholdt tem a limitação de não levar em consideração fatores como crises financeiras, por exemplo, de modo que o preço do ativo é reflexo exclusivamente da opinião dos agentes. Caso o mercado financeiro real fosse comportado desta forma, bastaria que todos os investidores tomassem a mesma opinião para que não houvesse prejuízo.

Porém, na vida real, a imprevisibilidade de ativos acaba compensando alguns agentes a terem opinião contrária à maioria, de modo que, quando um ativo em alta cai na bolsa, estes ganham um retorno consideravelmente maior do que os que optaram por partilhar da mesma opinião.

Na figura abaixo, têm-se uma malha que representa o sistema, onde cada pixel representa um agente, e as cores indicam qual a opinião de cada agente ao longo do tempo. Foram utilizadas as mesmas condições iniciais que no exemplo acima.

Comportamento dos agentes ao longo do tempo, com o mesmo número de agentes para cada opinião em t=0.	Comportamento dos agentes ao longo do tempo, com 100% de agentes de opinião 3 em t=0.	Comportamento dos agentes ao longo do tempo, com 100% de agentes de opinião 1 em t=0.	Comportamento dos agentes ao longo do tempo, com 100% de agentes de opinião 2 em t=0.

A figura reforça o que já foi confirmado antes, sobre o sistema tender a colapsar inteiramente para a opinião 2. Mais um fator interessante que pode ser observado está na concentração de agentes de opinião 1 ao longo do tempo. Há um momento logo no início das simulações onde há a formação de pequenos clusters de opinião 1, que logo se dissolvem e viram de opinião 2.

No início da simulação, ainda há um breve retorno para agentes que apostam contra a maioria, o que justifica a formação destes aglomerados. Porém, com o passar do tempo e a maioria da rede sendo formada pela opinião 2, estes clusters somem, e a malha fica 100% coberta de opinião 2.

Conclusões

Tendo em vista o que foi abordado até aqui, algumas conclusões podem ser tiradas:

simulações com um número de agentes da ordem de 1024 (grade 32x32) descrevem melhor um mercado financeiro, visto que apresentam clusters de volatilidade e possuem uma distribuição de retornos próxima a de um índice de mercado;
através de simulações, verificou-se o perfil de investidores em mercados reais, baseado em suas opiniões. Em geral, a população que investe em ações é distribuída de forma que a grande maioria (aproximadamente 85%) atua de forma irracional (sem uma estratégia), uma parcela que gira em torno de 10% atua como fundamentalista (atua de acordo com a minoria do sistema) e aproximadamente 5% atua de maneira chartista (segue a opinião da maioria).
ao permitir que agentes troquem de opinião livremente durante a simulação, observou-se a tendência do sistema de colapsar para um estado em que todos ou quase todos os agentes sigam o que a maioria do sistema está fazendo. Esta situação não descreve bem um sistema real, já que em um sistema financeiro existem fatores imprevisíveis que compõem o preço de ativos, não somente a opinião dos agentes. A partir do momento em que há possibilidade da minoria lucrar, surgem agentes de opinião 1. Também observou-se que, ao "forçar" que um agente que perde troque de opinião, elimina-se do sistema a opinião 3, ou seja, todos os agentes tornam-se minimamente racionais.

Programas

Código para obter a série temporal do preço dos índices financeiros

import pandas_datareader as pdr
from datetime import datetime

ibov = pdr.get_data_yahoo(symbols='^BVSP',start=datetime(1995,1,1),end=datetime(2021,1,1))
SP500 =  pdr.get_data_yahoo(symbols='^GSPC',start=datetime(1995,1,1),end=datetime(2021,1,1))
DJ = pdr.get_data_yahoo(symbols='^DJI',start=datetime(1995,1,1),end=datetime(2021,1,1))
Nasdaq = pdr.get_data_yahoo(symbols='^IXIC',start=datetime(1995,1,1),end=datetime(2021,1,1))

Código para calcular os retornos e para os normalizar

def ret(x):
    N = len(x)

    y = []
    for i in range(N-1):
        r = np.log(x[i+1])-np.log(x[i])
        y.append(r)

    return y

def normalize(x):
    
    N = len(x)    
    y = []
    
    for i in range(N):
        n = x[i] - np.mean(x)
        n = n/np.std(x)
        
        y.append(n)
        
    return np.array(y)

Código que realiza a estimativa de densidade kernel

def gauss(x,mean,std_dev):
    u = (x - mean) / std_dev
    c = 1 / (np.sqrt(2 * np.pi))
    return c * np.exp(- 0.5 * u ** 2)

def kde(x,kernel="gauss",bw=0.1,n_points=1500):
    kernel_options = ["gauss"]
    data = np.array(x)
    x_kde = np.linspace(np.min(data)-bw,np.max(data)+bw,n_points)
    n = data.shape[0]     #Number of rows
    m = x_kde.shape[0]    #Number of columns
    kde_i = []
    if kernel == kernel_options[0]:
        for x in data:
            kde_i.append(gauss(x_kde,x,bw))
    else:
        print("Kernel not found!")
        print("Kernel options are:")
        for k in kernel_options:
            print(" - " + k)
        return np.nan
          
    kde_i = np.array(kde_i).reshape(n,m)
    
    kde = np.array([np.sum(kde_i[:,i]) for i in np.arange(m)])
    kde_norm = kde / np.sum(kde)
    return x_kde,kde_norm

Código para calcular a auto correlação das volatilidades

from scipy.signal import get_window
from scipy.fft import rfft, rfftfreq, irfft

def fft_calculation(y,window="parzen"):
    w = get_window(window=window, Nx=len(y))
    f = y * w
    freqs = np.fft.rfftfreq(len(y))
    return freqs,np.abs(rfft(f))

def autocorr_calculation(y):
    ift = irfft(np.abs(y) ** 2)
    ift_norm = np.abs(ift) / np.abs(ift).max()
    return ift_norm

w_ibov, s_ibov = fft_calculation(np.array(normalize(r_ibov))**2)
c_ibov = autocorr_calculation(s_ibov)

w_sp500, s_sp500 = fft_calculation(np.array(normalize(r_SP500))**2)
c_sp500 = autocorr_calculation(s_sp500)

w_DJ, s_DJ = fft_calculation(np.array(normalize(r_DJ))**2)
c_DJ = autocorr_calculation(s_DJ)

w_nasdaq, s_nasdaq = fft_calculation(np.array(normalize(r_nasdaq))**2)
c_nasdaq = autocorr_calculation(s_nasdaq)

w100, s100 = fft_calculation(np.array(normalize(R100))**2)
c100 = autocorr_calculation(s100)

w50, s50 = fft_calculation(np.array(normalize(R50))**2)
c50 = autocorr_calculation(s50)

w32, s32 = fft_calculation(np.array(normalize(R32))**2)
c32 = autocorr_calculation(s32)

w16, s16 = fft_calculation(np.array(normalize(R16))**2)
c16 = autocorr_calculation(s16)

Código para gerar a figura dos retornos

fig, ax = plt.subplots(4,2,figsize=(10,7))

ax[0][0].plot(np.arange(len(R16)),R16,'lightsteelblue',label='# agentes = 256')
ax[1][0].plot(np.arange(len(R32)),R32,'cornflowerblue',label='# agentes = 1024')
ax[2][0].plot(np.arange(len(R50)),R50,'blue',label='# agentes = 2500')
ax[3][0].plot(np.arange(len(R100)),R100,'midnightblue',label='# agentes = 10000')

ax[0][1].plot(np.arange(len(r_ibov)),normalize(r_ibov),'pink',label='Ibovespa')
ax[1][1].plot(np.arange(len(r_SP500)),normalize(r_SP500),'palevioletred',label='S&P500')
ax[2][1].plot(np.arange(len(r_DJ)),normalize(r_DJ),'mediumvioletred',label='Dow Jones')
ax[3][1].plot(np.arange(len(r_nasdaq)),normalize(r_nasdaq),'purple',label='NASDAQ')

# Setting labels & titles

fig.suptitle('Retornos(t)',fontsize=14)
fig.text(0.5,0, 't', ha='center', fontsize=12)
fig.text(0.28,0.93, 'Bornholdt', ha='center', fontsize=12)
fig.text(0.77,0.93, 'Índices', ha='center', fontsize=12)
fig.text(0, 0.5, 'r(t)', va='center', rotation='vertical',fontsize=12)


for aa in ax: 
    for a in aa:
        a.xaxis.set_major_locator(plt.MaxNLocator(8))
        a.yaxis.set_major_locator(plt.MaxNLocator(5))
        a.legend(loc='upper left')

fig.tight_layout()
plt.show()

Código para gerar a figura da distribuição dos retornos

y100 = np.array(R100)
y50 = np.array(R50)
y32 = np.array(R32)
y16 = np.array(R16)

kde100 = kde(y100,"gauss",np.std(y100)/3,1500)
kde50 = kde(y50,"gauss",np.std(y50)/3,1500)
kde32 = kde(y32,"gauss",np.std(y32)/3,1500)
kde16 = kde(y16,"gauss",np.std(y16)/3,1500)

y_ibov = np.array(normalize(r_ibov))
y_sp500 = np.array(normalize(r_SP500))
y_dj = np.array(normalize(r_DJ))
y_nasdaq = np.array(normalize(r_nasdaq))

kde_ibov = kde(y_ibov,"gauss",np.std(y_ibov)/3,1500)
kde_sp500= kde(y_sp500,"gauss",np.std(y_sp500)/3,1500)
kde_dj = kde(y_dj,"gauss",np.std(y_dj)/3,1500)
kde_nasdaq = kde(y_nasdaq,"gauss",np.std(y_nasdaq)/3,1500)

fig, ax = plt.subplots(figsize=(10,6))


ax.set_xlabel('Retorno normalizado',fontsize=12)
ax.set_ylabel('Densidade de probabilidade',fontsize=12)
fig.suptitle('Distribuição do retorno (normalizado)',fontsize=14)

plt.plot(kde100[0],kde100[1],label='# agentes 10000')
plt.plot(kde50[0],kde50[1],label='# agentes 2500')
plt.plot(kde32[0],kde32[1],label='# agentes 1024')
plt.plot(kde16[0],kde16[1],label='16')

plt.plot(kde_ibov[0],kde_ibov[1],'--',label='Ibovespa')
plt.plot(kde_sp500[0],kde_sp500[1],'--', label='S&P500')
plt.plot(kde_dj[0],kde_dj[1],'--', label='Dow Jones')
plt.plot(kde_nasdaq[0],kde_nasdaq[1],'--', label='NASDAQ')


ax.set_xlim(-3,3)
plt.legend()
fig.tight_layout()

Código para gerar a figura da auto correlação das volatilidades

fig, ax = plt.subplots(4,2,figsize=(10,7))

ax[0][0].plot(c16[:int(len(c16)/2)],'lightsteelblue',label='# agentes = 256')
ax[1][0].plot(c32[:int(len(c32)/2)],'cornflowerblue',label='# agentes = 1024')
ax[2][0].plot(c50[:int(len(c50)/2)],'blue',label='# agentes = 2500')
ax[3][0].plot(c100[:int(len(c100)/2)],'midnightblue',label='# agentes = 10000')

ax[0][1].plot(c_ibov[:int(len(c_ibov)/2)],'pink',label='Ibovespa')
ax[1][1].plot(c_sp500[:int(len(c_sp500)/2)],'palevioletred',label='S&P500')
ax[2][1].plot(c_DJ[:int(len(c_DJ)/2)],'mediumvioletred',label='Dow Jones')
ax[3][1].plot(c_nasdaq[:int(len(c_nasdaq)/2)],'purple',label='NASDAQ')

# Setting labels & titles
#ax[2].set_xlabel('Data',fontsize=12)
#ax[1].set_ylabel('Retorno',fontsize=12)
fig.suptitle('Autocorrelação das volatilidades',fontsize=14)
fig.text(0.5,0, '$\\tau$', ha='center', fontsize=12)
fig.text(0.28,0.93, 'Bornholdt', ha='center', fontsize=12)
fig.text(0.77,0.93, 'Indíces', ha='center', fontsize=12)
fig.text(0, 0.5, '$A(\\tau)$', va='center', rotation='vertical',fontsize=12)



for aa in ax: 
    for a in aa:
        a.xaxis.set_major_locator(plt.MaxNLocator(8))
        a.yaxis.set_major_locator(plt.MaxNLocator(5))
        a.legend(loc='lower left')
        a.set_yscale('log')
        a.set_xscale('log')
        

fig.tight_layout()
plt.show()

Objeto que implementa o Modelo de Bornholdt em Python

class Bornholdt():
    def __init__(self,N_rows,N_cols,T):
        self.L = N_rows     #Size of mesh in X
        self.N = N_cols     #Size of mesh in Y
        self.beta = 1/T     #beta ~ 1/T
        self.state = init_state(self.L,self.N)     #Initialize the network
        self.Mag = calc_mag(self.state)            #Magnetization for a specific state

    def update_spins(self):
        """ Provides a single update in the network """
        self.state = heat_bath(self.state,self.beta)
        self.Mag = calc_mag(self.state)

    def reset_state(self):
        """ Reset the grid to the initial state """
        self.state = init_state(self.L,self.N)
        self.Mag = calc_mag(self.state)
        
    def magnetization(self,mcSteps,plot=False,number_prints=8):
        """ Calculates the magnetization and updates the spins of system. """
        n = int(mcSteps/number_prints)
        Magnetism = np.zeros(mcSteps) 
        #Heat-Bath
        for j in np.arange(mcSteps):
            if j % n == 0:
                print("MC Sweep nº: {}".format(j))
            self.update_spins()
            Magnetism[j] = self.Mag
        #Plot
        if plot:
            fig,ax = plt.subplots(1,1,figsize=(8,4))
            ax.set_ylabel("M(t)",fontsize=12)
            ax.set_xlabel("t",fontsize=12)
            ax.plot(Magnetism,lw=0.8,color='black')
            ax.set_xlim(0,len(Magnetism))
            plt.tight_layout()
            plt.show()
        return Magnetism
    
    def print_state(self):
        """ Plots the current status of the grid, in a 2D mesh """
        X,Y = np.meshgrid(np.arange(self.L),np.arange(self.N))
        
        fig,ax = plt.subplots(1,1,figsize=(3,3))
        ax.set_xticks([])
        ax.set_yticks([])
        ax.pcolor(X, Y, self.state, cmap=plt.cm.Greys,shading='auto',linewidth=0,rasterized=True)
        
        plt.tight_layout()
        plt.show()
    
    def plot_grid(self,mcSteps,n_plots,number_prints=8):
        """ Built a grid with 2D meshs plots, perfect to visualize the system. """
        n = n_plots**2
        X,Y = np.meshgrid(np.arange(self.L),np.arange(self.N))
        fig,ax = plt.subplots(n_plots,n_plots,figsize=(8,8))
        print_graph = int(mcSteps/n)
        snapshots = []
        #Constructing snapshots list
        for i in np.arange(mcSteps):
            if (i % print_graph==0) and (len(snapshots) < n):
                snapshots.append(i)
        #Heat-Bath
        count = 0
        for j in np.arange(mcSteps):
            if j % int(mcSteps/number_prints) == 0:
                print("-- MC sweep nº {}...".format(j))
            if j in snapshots:
                #Plot config.
                x_i = int(count/n_plots)%n_plots
                y_i = int(count%n_plots)
                ax[x_i,y_i].pcolor(X, Y, self.state, cmap=plt.cm.Greys,shading='auto',linewidth=0,rasterized=True)
                ax[x_i,y_i].set_title("t = {}".format(j),fontsize=22)
                ax[x_i,y_i].set_xticks([])
                ax[x_i,y_i].set_yticks([])
                count += 1
            self.update_spins()
        #Plot config.
        plt.tight_layout()
        plt.show()

Implementação do Modelo de Bornholdt com variação de opiniões ao longo da simulação (Julia)

using Statistics
println("Statistics imported sucessfully")

function calc_ret(M)
    N = length(M)
    corrected_M = []
    for m in M
        if m == 0.0
            append!(corrected_M,0.000001)
        else
            append!(corrected_M,m)
        end
    end
    r = []
    for i in 2:N
        ret = log(abs(corrected_M[i]/corrected_M[i-1]))
        append!(r,ret)
    end
    r_ = mean(r)
    sigma_r = std(r)
    r_norm = [(ret-r_)/sigma_r for ret in r]
    return r_norm
end

function agent_opinion(p1,p2)
    opinion = rand()
    if opinion <= p1
        return +1,+1.0 #Strategy 1 - Player that wants to join the minority (fundamentalist)
    elseif opinion <= (p2+p1)
        return +2,-1.0 #Strategy 2 - Player wants to join majority of agents (chartist)
    else
        return +3,(2*rand(0:1) - 1) #Strategy 3 - Dumb irrational player
    end
end

function flip_opinion(current_opinion)
    opinions = [1,2,3]
    d = Dict()
    d[1] = +1.0
    d[2] = -1.0
    d[3] = (rand() - 0.5)
    splice!(opinions,Int64(current_opinion))
    n = length(opinions)
    new_opinion = rand(1:n)
    return new_opinion,d[new_opinion]
end

function count_states(opinion_mesh)
    N,L = size(opinion_mesh) 
    opinion1 = 0
    opinion2 = 0
    opinion3 = 0
    for i in 1:N
        for j in 1:L
            if opinion_mesh[i,j] == 1
                opinion1 += 1/(N*L)
            elseif opinion_mesh[i,j] == 2
                opinion2 += 1/(N*L)
            else
                opinion3 += 1/(N*L)
            end
        end
    end
    return [opinion1,opinion2,opinion3]
end

function get_C(spins,p1,p2)
    N,L = size(spins)
    C = zeros(N,L)
    opinions = [0,0,0]
    opinion_mesh = zeros(N,L)
    for i in 1:N
        for j in 1:L
            opinion, C_value = agent_opinion(p1,p2)
            C[i,j] = C_value
            opinion_mesh[i,j] = opinion
            opinions[opinion] += 1
        end
    end
    return opinions,C,opinion_mesh
end

function init_network(N,L)
    spins = zeros(N,L)
    for i in 1:N
        for j in 1:L
            r = rand()
            if r < 0.5
                spins[i,j] = -1
            else
                spins[i,j] = +1
            end
        end
    end
    return spins
end

function find_nb(spins,a,b)
    N,L = size(spins)
    #Adjusting for a
    if (a-1) <= 0
        back_a = N
    else
        back_a = a-1
    end

    if (a+1) > N
        forward_a = 1
    else
        forward_a = a+1
    end

    #Adjusting for b
    if (b-1) <= 0
        back_b = N
    else
        back_b = b-1
    end

    if (b+1) > L
        forward_b = 1
    else
        forward_b = b+1
    end

    left = spins[a,back_b]
    right = spins[a,forward_b]
    top = spins[back_a,b]
    bottom = spins[forward_a,b]
    return [left,right,top,bottom]
end    

function calc_spin(arg)
    p = 1 / (1 + exp((-1)*arg))
    r = rand()
    if r < p
        spin = +1
    else
        spin = -1
    end
    return spin
end

function heat_bath(spins,beta,alpha,C,change_opinions)
    N,L = size(spins)
    J = 1
    a_vec = rand(1:N,N^2)
    b_vec = rand(1:L,L^2)
    for (a,b) in zip(a_vec,b_vec)
        s = spins[a,b]
        nb = sum(find_nb(spins,a,b))
        local_field_h = J*nb - alpha*C[a,b]*sum(spins)/(N*L)
        spins[a,b] = calc_spin(2*beta*local_field_h)
        if opinion_mesh[a,b] == 3
            C[a,b] = (2*rand(0:1) - 1)
        end    
        if change_opinions
            if s != spins[a,b]
                (opinion,C_value) = flip_opinion(opinion_mesh[a,b])
                C[a,b] = C_value #If the agent changes signal, he changes strategy
                opinion_mesh[a,b] = opinion
            end
        end
    end
    return spins
end

function calc_mag(spins)
    N,L = size(spins)
    m = sum(spins)
    return m / (N*L)
end

function write_to_file(filename,data1,data2)
    N = length(data1)
    open(filename,"w") do file
        write(file,"M(t),r(t)\n")
        write(file,(string(data1[1])*",-\n"))
        for i in 2:N
            write(file,string(data1[i])*","*string(data2[i-1])*"\n")
        end
    end
end

function print_state(filename,C)
    N,L = size(C)
    file = open(filename,"w")
    for i in 1:N
        for j in 1:L
            write(file,string(C[i,j])*",")
        end
        write(file,"\n")
    end
    close(file)
end

function print_opinions(filename,opinions)
    N,L = size(opinions)
    file = open(filename,"w")
    for i in 1:N
        for j in 1:L
            write(file,string(opinions[i,j])*",")
        end
        write(file,"\n")
    end
end

N = 32
L = 32
mcSteps = 3000
n_prints = Int(floor(mcSteps/4))
T = 1.5
beta = 1/T
alpha = 4.0

p1 = 1/3 # Prob. opinion #1: fundamentalist
p2 = 1/3 # Prob. opinion #2: chartist
change_opinions = true # If true, opinion of agents can change during simulation.

spins = init_network(N,N)
M = zeros(mcSteps)
opinions,C,opinion_mesh = get_C(spins,p1,p2)

println("Network opinions:")
for (i,opinion) in enumerate(opinions)
    println("Opinion ",i,": ",100*opinion/(N*L)," %")
end
println("====================================")
perc_opinions = zeros(mcSteps,3)
plot_grid = floor(mcSteps/15)
shots = [1]
append!(shots,[i*plot_grid for i in 1:15])
count = 0

for i in 1:mcSteps
    if (i in shots)
        println("----- Grid in ",i)
        filename = "state_MC"*string(count)*".txt"
        print_state(filename,opinion_mesh)
        global count += 1
    end
    percs = count_states(opinion_mesh)
    for j in 1:3
        perc_opinions[i,j] = percs[j]
    end
    heat_bath(spins,beta,alpha,C,change_opinions)
    M[i] = calc_mag(spins)
    if i % n_prints == 0
        println("-- Sweep nº ",i," concluded.")
    end
end
filename = "mag_ret_N"*string(N*L)*".txt"
write_to_file(filename,M,calc_ret(M))
print_opinions("opinions_N"*string(N^2)*".txt",perc_opinions)

Implementação do Modelo de Bornholdt simples em C

OBS: código não otimizado.

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>

#define N 32
#define L 32
#define MAX_STRING_SIZE 15
#define NUMBER_OF_PLOTS 9
int spins[N][L];

double rng(){
    return (float)rand()/RAND_MAX;
}

int int_rng(int max_value){
    return rand() % max_value;
}

void initialize_network(){
    for(int i=0;i<N;i++){
        for(int j=0;j<L;j++){
            spins[i][j]=2*int_rng(2)-1;
        }
    }
}

double sum_array(){
    double sum = 0.0;
    for(int i=0;i<N;i++){
        for(int j=0;j<L;j++){
            sum = sum + spins[i][j];
        }
    }
    return sum;
}

int calc_spin(double arg){
    double p = 1 / (1 + exp(arg));
    double r = rng();
    if(r < p){
        return +1;
    }
    else{
        return -1;
    }
}

double calc_magnetization(){
    double m = 0.0;
    for(int i=0;i<N;i++){
        for(int j=0;j<L;j++){
            m = m + spins[i][j];
        }
    }
    return (double) (m/(L*N));
}

void print_state(char filename[]){
    FILE *grid;
    grid = fopen(filename,"w");
    for(int i=0;i<N;i++){
        for(int j=0;j<L;j++){
            fprintf(grid,"%d,",spins[i][j]);
        }
        fprintf(grid,"\n");
    }
    fclose(grid);
}

void heat_bath(double beta){
    double alpha = 4.0;
    double J = 1.0;
    int a;
    int b;
    int s;
    int nb_left,nb_right,nb_top,nb_bottom,sum_nb;
    double arg;
    for(int i=0;i<N;i++){
        for(int j=0;j<L;j++){
            a = int_rng(N);
            b = int_rng(L);
            s = spins[a][b];
            nb_top = spins[(a+N-1)%N][b];
            nb_bottom = spins[(a+1)%N][b];
            nb_left = spins[a][(b+L-1)%L];
            nb_right = spins[a][(b+1)%L];
            sum_nb = nb_left+nb_top+nb_right+nb_bottom;
            arg = (double)2*beta*(J*sum_nb - alpha*s*abs(sum_array(*spins))/(L*N));
            spins[a][b] = calc_spin(-arg);
        }
    }
}

void main(){
    srand(time(NULL));
    //!Declarate variables
    clock_t begin;
    clock_t end;
    double T = 1.5;
    double beta = 1/T;
    int mcSteps = 1000;
    int n_prints = (int)(mcSteps/5);
    char filename[][MAX_STRING_SIZE] = {"grid1.txt","grid2.txt","grid3.txt","grid4.txt","grid5.txt","grid6.txt","grid7.txt","grid8.txt","grid9.txt"};
    int snapshots[NUMBER_OF_PLOTS];
    for(int i=0;i<NUMBER_OF_PLOTS;i++){
        snapshots[i] = i*(mcSteps/NUMBER_OF_PLOTS);
    }
    int plot_index = 0;
    FILE *mag_file;
    mag_file = fopen("magnetization.dat","w");
    //////////////////////////////////////////////////////////////////////!
    begin = clock();
    //!Init. Network
    initialize_network();
    end = clock();
    printf("Network initiated\n");
    printf("Executed time: %lf s\n", (double)(end-begin)/CLOCKS_PER_SEC);
    printf("==========================================================\n");
    //////////////////////////////////////////////////////////////////////!
    begin = clock();
    for(int j=0;j<mcSteps;j++){
        if(j % n_prints == 0){
            printf("Simulation %d...\n",j);
        }
        if(j == snapshots[plot_index]){
            printf("-- Building %s...\n",filename[plot_index]);
            print_state(filename[plot_index]);
            plot_index++;
        }
        //!Heat-Bath
        heat_bath(beta);
        //!Magnetization
        fprintf(mag_file,"%lf\n",calc_magnetization());
    }
    fclose(mag_file);
    end = clock();
    if(N>40){
        printf("Executed time: %lf min\n", (double)(end-begin)/(60*CLOCKS_PER_SEC));
    }
    else{
        printf("Executed time: %lf s\n", (double)(end-begin)/CLOCKS_PER_SEC);
    }
    printf("==========================================================\n");
    system("python plot.py");
}

Referências

↑ ^1,0 ^1,1 Bornholdt, Stefan. (2011). Expectation bubbles in a spin model of markets: Intermittency from frustration across scales. International Journal of Modern Physics C. 12. 10.1142/S0129183101001845.
↑ https://en.wikipedia.org/wiki/Ising_model
↑ Lux, Thomas & Marchesi, Michele. (1998). Scaling and Criticality in a Stochastic Multi-Agent Model of a Financial Market. Nature. 397. 10.1038/17290.
↑ ^4,0 ^4,1 Lux, Thomas & Marchesi, Michele. (1998). Volatility Clustering in Financial Markets: A MicroSimulation of Interacting Agents. International Journal of Theoretical and Applied Finance. 3. 10.1142/S0219024900000826.
↑ Retornos. Portal Action.
↑ ^6,0 ^6,1 Kernel Density Estimation. Wikipedia
↑ Da Cunha, Carlo & Silva, Roberto. (2019). Relevant Stylized Facts About Bitcoin: Fluctuations, First Return Probability, and Natural Phenomena.
↑ Wiener-Khinchin Theorem. Wolfram Math World
↑ Wiener–Khinchin theorem. Wikipedia

[modelo-1] 1,0 ^1,1 Bornholdt, Stefan. (2011). Expectation bubbles in a spin model of markets: Intermittency from frustration across scales. International Journal of Modern Physics C. 12. 10.1142/S0129183101001845.

[ising-2] ttps://en.wikipedia.org/wiki/Ising_model

[bonrholdt12-3] Lux, Thomas & Marchesi, Michele. (1998). Scaling and Criticality in a Stochastic Multi-Agent Model of a Financial Market. Nature. 397. 10.1038/17290.

[bonrholdt13-4] 4,0 ^4,1 Lux, Thomas & Marchesi, Michele. (1998). Volatility Clustering in Financial Markets: A MicroSimulation of Interacting Agents. International Journal of Theoretical and Applied Finance. 3. 10.1142/S0219024900000826.

[ret-5] Retornos. Portal Action.

[kde-6] 6,0 ^6,1 Kernel Density Estimation. Wikipedia

[fatoestilizado-7] Da Cunha, Carlo & Silva, Roberto. (2019). Relevant Stylized Facts About Bitcoin: Fluctuations, First Return Probability, and Natural Phenomena.

[wiener1-8] Wiener-Khinchin Theorem. Wolfram Math World

[wiener2-9] Wiener–Khinchin theorem. Wikipedia

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9]