Taller: Modelado de agentes económicos con Python.

Seminario de Economía y Complejidad (CEIICH-UNAM), 8 de abril, 2025

Actividad: Implementación de Modelo de Agentes de Dragulescu A. & Yakovenko V. (2000)[1].

Instructor: Dr. Andrés García Medina

email: andgarm.n@gmail.com

sitio: https://sites.google.com/view/andresgm/home.

[1] Dragulescu, A., & Yakovenko, V. M. (2000). Statistical mechanics of money. The European Physical Journal B-Condensed Matter and Complex Systems, 17, 723-729.

import numpy as np
import matplotlib.pyplot as plt

# Parametros
N =  1000 # Número de agentes
M = 10000  # Dinero total en el sistema
iterations = 500000  # Número de transacciones

# Inicializar la distribución de dinero: distribución equitativa del dinero total entre los agentes
money = np.ones(N) * (M / N)

# Proceso de intercambio de dinero (aleatoria)
for _ in range(iterations):
    # Elegir dos agentes aleatorios
    i, j = np.random.randint(0, N, 2)

    # Asegurar de que no son el mismo agente
    if i != j:
        # Fracción aleatoria de dinero a intercambiar
        epsilon = np.random.random()

        # Calcular la cantidad de dinero a transferir de i a j
        delta_m = epsilon * (M/N)

        # Asegúrese de que el dinero se mantenga no negativo
        if money[i] >= delta_m:
            money[i] -= delta_m
            money[j] += delta_m

# Grafique el histograma de la distribución del dinero
plt.hist(money, bins=50, density=True, alpha=0.6, color='g', label='Datos simulados')

#Ajuste a la distribución de Boltzmann-Gibbs (ajuste exponencial)
average_money = M / N
x = np.linspace(0, max(money), 1000)
boltzmann_gibbs = (1 / average_money) * np.exp(-x / average_money)

plt.plot(x, boltzmann_gibbs, 'r-', linewidth=2, label='Distribución Boltzmann-Gibbs')

plt.xlabel('Dinero por Agente')
plt.ylabel('Densidad de probabilidad')
plt.title('Distribución de dinero entre agentes')
plt.legend()
plt.grid(True)
plt.savefig('Yakovenko_simulation.png')
plt.show()

Podemos integrar más reglas de intercambio de dinero

# Función para seleccionar un par de agentes aleatorios y realizar el intercambio
def exchange_money(money, rule_type=3):
     # Elegir dos agentes aleatorios
    i, j = np.random.randint(0, N, 2)

    # Asegurar de que no son el mismo agente
    if i != j:


        # Regla 1: intercambio de una cantidad fija de dinero
        if rule_type == 1:
            delta_m = 1
            if money[i] >= delta_m:
                money[i] -= delta_m
                money[j] += delta_m
        # Regla 2: intercambio de una fracción aleatoria del dinero promedio
        elif rule_type == 2:
            nu = np.random.random()
            delta_m = nu * (money[i] + money[j]) / 2
            if money[i] >= delta_m:
                money[j] += delta_m
                money[i] -= delta_m
        # Regla 3: intercambio de una fracción aleatoria del dinero promedio en el sistema
        elif rule_type == 3:
            nu = np.random.random()
            delta_m = nu * M / N
            if money[i] >= delta_m:
                money[i] -= delta_m
                money[j] += delta_m
    return(money)

Mostramos la distribución de la riqueza cada 100,000 pasos para alguna de las reglas de intercambio

# Parametros
N =  1000 # Número de agentes
M = 10000  # Dinero total en el sistema
iterations = 500001  # Número de transacciones

# Inicializar la distribución de dinero: distribución equitativa del dinero total entre los agentes
money = np.ones(N) * (M / N)

# Simulación
for step in range(iterations):
    money = exchange_money(money, rule_type=3)

    # Cada 100000 pasos, se puede visualizar la distribución
    if step % 100000 == 0:
        plt.clf()
        plt.close()
        plt.figure()
        # Grafique el histograma de la distribución del dinero
        plt.hist(money, bins=50, density=True, alpha=0.6, color='g', label='Datos simulados')

        #Ajuste a la distribución de Boltzmann-Gibbs (ajuste exponencial)
        average_money = M / N
        x = np.linspace(0, max(money), 1000)
        boltzmann_gibbs = (1 / average_money) * np.exp(-x / average_money)

        plt.plot(x, boltzmann_gibbs, 'r-', linewidth=2, label='Distribución Boltzmann-Gibbs')

        plt.xlabel('Dinero por Agente')
        plt.title(f'Distribución de Dinero en el Paso {step}')
        plt.legend()
        plt.xlim(0, max(money))
        plt.grid(True)
        plt.pause(0.01)

# Mostrar la distribución final
plt.show()