Primeira vez aqui? Seja bem vindo e cheque o FAQ!
x

Explorando a base de dados "Strategy Choice in The Infinitely Repeated Prisoners’ Dilemma" de Pedro Dal Bó and Guillaume R. Frechette (2019), AER

0 votos
17 visitas
perguntada Mai 18 em Economia por Ricardo Saldanha (1 ponto)  
editado Mai 21 por Ricardo Saldanha

Explorando a base de dados principal do artigo "Strategy Choice in The Infinitely Repeated Prisoners’ Dilemma" de Pedro Dal Bó and Guillaume R. Frechette - AMERICAN ECONOMIC REVIEW, VOL. 109, NO. 11, NOVEMBER 2019 (pp. 3929-52).
O artigo e as bases estão disponíveis em:
https://www.aeaweb.org/articles?id=10.1257/aer.20181480.
A base principal chama-se "datstratfinal_ac.dta". Há bases auxiliares, que servem para testes de robustez e para comparação com os resultados de um artigo anterior dos autores. Mas nos concentraremos na principal.

Compartilhe

1 Resposta

0 votos
respondida Mai 18 por Ricardo Saldanha (1 ponto)  
editado Mai 19 por Ricardo Saldanha

Objetivo principal do artigo
Os autores estão interessados em estimar prevalências de estratégias que surgem em jogos repetidos do Dilema do Prisioneiro. Para isso pedem aos participantes que elicitem estratégias; i.e., tracem um plano.

Motivação de explorar a base
O estudo de jogos repetidos tem ampla aplicação, como em organização industrial e contratos informais.
Trata-se de um experimento randomizado, de forma que é possível utilizar um método mais direto, que é a Estimação por Máxima Verossimilhança (MLE), sem perda de rigor.

Base econômica
Estimação de prevalências de estratégias. Os resultados podem ser comparados com o previsto pela teoria e direciona estudos para as estratégias mais relevantes.
Há dificuldades associadas a estimar estratégias em jogos infinitos apenas observando o comportamento: conjunto infinito de estratégias possíveis; não se observa o que os indivíduos teriam feito em outro cenário.
A abordagem de MLE já havia sido usada (Dal Bó and Fréchette 2011; Fudenberg, Rand, and Dreber 2012). Mas não com elicitação.
Embora a elicitação já tenha sido examinada, em geral são casos finitos ou o foco não é a prevalência de estratégias. Romero and Rosokha (2018) são uma exceção; porém, o método de elicitação é mais livre, podendo surgir problemas como planos inconsistentes.

Base de dados - experimento
O jogo, com probabilidade \(\delta\) de continuar:
\begin{array}{lcc}
& C & D \\ \hline
C & R,R & 12,50 \\ \hline
D & 50,12 & 25,25 \\
\end{array}
R será 32 ou 48 e \(\delta\) será 1/2, 3/4 ou 9/10, conforme o tratamento.
São três fases:
Fase 1: Várias partidas (superjogos) entre pares sorteados, sem plano de ação (acaba com o fim da primeira partida encerrada após 20 min).
Fase 2: Como a Fase 1, mas os participantes devem traçar um plano de ação e atualizá-lo a cada rodada (exatamente 20 min, podendo interromper uma partida no meio).
Fase 3: Escolhas são interrompidas e o computador joga sozinho mais 15 rodadas conforme a última estratégia definida.
Abaixo o quadro que aparece para os participantes. Os indivíduos escolhem ações para cada uma das 4 possibilidades do jogo.
A imagem será apresentada aqui.
Em algumas seções, os autores apresentam um menu de estratégias (indicado na base com a binária "menu") antes da opção de montar estratégia personalizada como acima. Não se encontram, porém, significativas diferenças entre uma forma ou outra de elicitar.
As escolhas mais presentes na amostra são compatíveis com as estratégias simples AC, AD, GRIM, TFT, STFT (TFT começando sem cooperar), WSLS ("win stay, loose shift"). E são essas estreatégias que farão parte do modelo, para extrapolar a amostra.

Modelo a ser reproduzido
\( {\displaystyle Pr_{i}(s^{k}) = \prod_{M_{i}}\prod_{R_{im}}(\beta)^{I_{imr}^{k}}(1-\beta)^{1-I_{imr}^{k}} = \beta^{Y_{i}}(1-\beta)^{N_{i}-Y_{i}} }\)
\(I_{imr}^{k}\) é dummy=1 se o participante i jogou como prevê a estratégia k na rodada r (da partida m). $Yi$ é o número de sucessos (soma de \(I_i\)) entre o total \(N_{i}\) de rodadas jogadas pelo indivíduo.
Estimação por Máxima Verossimilhança:
\( (\hat{\beta},\hat{\phi^{1}},...,\hat{\phi^{K}}) \) é tal que satisfaz
\( {\displaystyle \max
{\beta \in (\frac{1}{2},1),\phi\geq0} \sum{I}ln\bigg{(}\sum{K}\phi^{k}Pr_{i}(s^{k})\bigg{)} \\
\text{s.t. } \phi^{1}+...+\phi^{K}=1 } \)

Para entender a interpretação dos parâmetros:
Um \(\beta=0.8\) significa que quando o plano declarado pelo indivíduo é compatível com uma dada estratégia, ele jogará de acordo com essa estratégia 80% das vezes.
Tome \(\phi^{AD}=0.6\), \(\phi^{GRIM}=0.3\) e \(\phi^{TFT}=0.1\). Então modelo entende que em 60% partidas é adotado AD, 30% GRIM e 10% TFT.
Para o 1º round, e.g., \(s^{AD}=D\) e \(s^{GRIM}=s^{TFT}=C\).
Por fim, a taxa de cooperação esperada no 1º round é \(0.8*0.4+(1-0.8)*(1-0.4)=0.44\).

Resultados
Apresentarei aqui apenas a tabela principal replicada, que é a Table 5 do artigo. Esta responde à pergunta na qual estamos interessados: se é possível criar um bom modelo de prevalência de estratégias. A rotina abaixo gera algumas outras tabelas, conforme explicado nos comentários.
A imagem será apresentada aqui.
Os resultados das prevalências e coeficientes são bem próximos dos apresentados pelos autores, com exceção talvez das primeiras estratégias do quarto tratamento.
As variância, por outra lado, ficaram em geral menores que a dos autores, identificando mais coeficientes estatisticamente diferentes do que no artigo. Ainda assim, os resultados não são tão diferentes.
A variância dos coeficientes é estimada por um bootstrap de 1000 iterações para cada tratamento; cada iteração reestima os coeficientes. Para obter uma subamostra, os autores mantêm todos os indivíduos e sorteiam partidas. Enquanto que eu sorteei subamostras de indivíduos. Isso provalmente explicação parte das diferenças nos resultados.

Código
Como explicado no arquivo "analyis\readme.txt" dos suplementes, parte do código é feito em Stata e parte em Matlab. A Seção IV, que tem a maior parcela do conteúdo replicado aqui, foi feita pelos autores em Matlab. Eu não consegui rodar o código do artigo no programa de código aberto Octave. Portanto, decidi elaborar meu próprio código, em Python, com base nas explicações do artigo.
São duas rotinas. A primeira faz manipulações básicas e cria tabelas de prevalências de estratégias na amostra. A segunda é responsável pelas estimações.

Tratamento prévio da base/prevalências amostrais:

# -*- coding: utf-8 -*-

import os
import pandas as pd
import numpy as np

os.chdir(r'filepath')

if __name__=="__main__":

    strats = ['ac','tft','wsls','grim','stft','ad']

    ''' [Ref. à Tabela 3 - prevalências]'''

    df = pd.read_stata('data/dat_strat_final_ac.dta') # Base original do experimento principal.
    df = df[df['phase']==2]
    cond1 = (df['last']==1) & (df['round']==1)
    df['strat_final'] = np.where(cond1,'other',np.nan) # np.where aqui faz a função de cond() do Stata.
    for strat in strats:
        df.loc[cond1 & (df[strat]==1), 'strat_final'] = strat

    # df[cond1] é o base com o filtro cond1. groupby()
    tab_prevalence = df[cond1].groupby('treatment')['strat_final'].value_counts(normalize=True)
    tab_prevalence = tab_prevalence.unstack().fillna(0)
    print(tab_prevalence)
    tab_prevalence.to_excel("Tab3_prevalence.xlsx")

    ''' [Ref. à Tabela 5 - prevalências] 
    Retiramos o primeiro quarto dos rounds, 
    considerado pelos autores como um período de aprendizado dos participantes.'''

    df = pd.read_stata('data/dat_strat_final_ac.dta')
    df = df[df['phase']==2]

    cond1 = df['round']==1
    df['strat_match'] = np.where(cond1,'other',np.nan) # np.where aqui faz a função de cond() do Stata.
    for strat in strats:
        df.loc[(cond1) & (df[strat]==1), 'strat_match'] = strat

    keep = ['treatment','date','id','match','round','coop','ocoop','strat_coop','menu','first','last','strat_match']
    df = df[keep]

    df = df.groupby('id').apply(lambda x: x.sort_values(['match', 'round'])).reset_index(drop=True)
    df['index_totrounds'] = df.groupby('id').cumcount() + 1
    df['num_totrounds'] = df.groupby('id')['round'].transform('count')
    # Eliminando o primeiro quarto das observações, cfe. comentado acima.
    filt = df['index_totrounds']>df['num_totrounds']/4
    df = df[filt]
    cond1 = df['round']==1
    tab_prevalence = df[cond1].groupby('treatment')['strat_match'].value_counts(normalize=True).unstack().fillna(0)   
    print(tab_prevalence)
    tab_prevalence.to_excel("Tab5_prevalence.xlsx")



    ''' [Ref. à Tabela 5 - base estimação]
     Queremos comparar o que as principais estratégias determinariam para dada jogada com o que o participando de fato jogou.'''

    # coop é dummy==1 para jogada de cooperação. Criando variável defasada, i.e., coop no round anterior da mesma partida.
    df['coop_m1'] = df.groupby(['id','match'])['coop'].shift(1)
    df['ocoop_m1'] = df.groupby(['id','match'])['ocoop'].shift(1) # Jogada defasada do adversário.

    # Variáveis com as jogadas previstas pelas estratégias, de acordo com a definição de cada estratégia.
    df['coop_ac'] = 1
    df['coop_ad'] = 0
    df['coop_grim'] = np.where(df['round']==1,1,np.where((df['coop_m1']==1) & (df['ocoop_m1']==1),1,0))
    df['coop_tft'] = np.where(df['round']==1,1,np.where(df['ocoop_m1']==1,1,0))
    df['coop_stft'] = np.where(df['round']==1,0,np.where(df['ocoop_m1']==1,1,0))
    df['coop_wsls'] = np.where(df['round']==1,1,np.where(df['coop_m1']==df['ocoop_m1'],1,0))

    # Dummies para jogada feita de acordo com a respectiva estratégia.
    newvars = ['success_'+strat for strat in strats]
    coopvars = ['coop_'+strat for strat in strats]
    for i in range(len(newvars)):
        df[newvars[i]] = np.where(df['coop']==df[coopvars[i]],1,0)

    # Antes de consolidar sucessos por participante, manteremos uma coluna que informa a escolha final do mesmo. 
    cond = (df['last']==1) & (df['round']==1) # round==1 para que fique só um registro por indivíduo.
    df['strat_final'] = np.where(cond,df['strat_match'],np.nan)

    df['num_trials'] = df.groupby('id')['round'].transform('count')
    # Base só com os dados filtrados por cond e colunas selecionadas; mesclado com as somas de successo, formará o consolidado. 
    df1 = df.loc[cond,['date','id','num_trials','strat_final','treatment','menu']]

    # Consolidando número de sucessos por participante; i.e., df2 é a tabela com as somas de sucesso por indivíduo.
    df2 = df.groupby('id')[newvars].sum()
    dic_newvars = {'success_'+strat:'num_success_'+strat for strat in strats}
    df2.rename(columns = dic_newvars, inplace = True)

    df = df1.merge(df2, on=['id'], how='inner') 
    df.to_csv('dataframe.csv', sep=';',index=False)

Estimações:

# -*- coding: utf-8 -*-

import os
import pandas as pd
import numpy as np
from scipy.optimize import minimize
from scipy.stats import chi2

os.chdir(r'filepath')
df_alltreats = pd.read_csv('dataframe.csv', sep=';') # gerado na rotina acima
df_alltreats.drop(columns=['date','strat_final'],inplace=True)

# Função de verossimilhança usando os parâmetros as init_guess em MLE().
def LL(params):
    beta = params[0]

    nums_success = ['num_success_'+x for x in strats]
    Pr = []
    for k in range(len(strats)):
        Pr_k = beta**df[nums_success[k]] * \
            (1-beta)**(df['num_trials']-df[nums_success[k]])
        Pr.append(Pr_k)

    loglik = -np.sum(np.log(np.dot(params[1:],Pr)))

    return loglik

def constraint(params):
    return np.sum(params[1:])-1


def MLE(strats,beta_guess,cte_phi_guess,bnd_beta,bnds_phi):

    init_guess = [beta_guess,]+[cte_phi_guess]*len(strats)
    bnds = (bnd_beta,) + (bnds_phi,)*len(strats)
    const = {'type':'eq', 'fun': constraint} # ref. à f. constraint acima

    results = minimize(LL,init_guess,method='SLSQP', bounds=bnds, \
                       constraints=const)

    dic_results = {}
    keys = ['beta']+[x for x in strats]
    for i in range(len(keys)):
        dic_results[keys[i]] = [results.x[i]]

    return dic_results


# Não funcionará com quaisquer dicionários; values têm que ser listas.
def merge_dics(dic_orig,dic_to_add):
    if bool(dic_orig): # Retorna False se o dicinário for vazio.
        keys_orig = list(dic_orig.keys())
        keys_to_add = list(dic_to_add.keys())
        for key in keys_to_add:
            if key in keys_orig:
                dic_orig[key] += dic_to_add[key]
            else:
                dic_orig[key] = dic_to_add[key]
    else:
        dic_orig = dic_to_add

    return dic_orig


def MLE_treats(cd_treatments):
    global df
    global strats

    # 1ª estimação
    dic_results = {}
    positive_coef_strats = []
    for treat in cd_treatments:
        filt_treat = (df_alltreats['treatment']==treat)
        df = df_alltreats[filt_treat].reset_index(drop=True)

        dic_results_row = MLE(strats,beta_guess,cte_phi_guess,bnd_beta,bnds_phi)
        dic_results = merge_dics(dic_results,dic_results_row)

        considerado_nulo=10**-8
        keys_positive = [key for key,value in dic_results_row.items() \
                                    if value[0]>considerado_nulo]
        positive_coef_strats.append(keys_positive[1:]) # [1:] p/ excluir o beta.

    tab1 = pd.DataFrame(dic_results)

    # 2ª estimação
    old_strats = strats
    dic_results = {}
    for treat in cd_treatments:
        filt_treat = (df_alltreats['treatment']==treat)
        df = df_alltreats[filt_treat].reset_index(drop=True)
        strats = positive_coef_strats[treat-1]
        zero_coef_strats = [x for x in old_strats if not x in strats]

        dic_results_row = MLE(strats,beta_guess,cte_phi_guess,bnd_beta,bnds_phi)
        # Adicionando colunas vazias, para o DataFrame.
        for x in zero_coef_strats:
            dic_results_row[x] = ['dropped']
        dic_results = merge_dics(dic_results,dic_results_row)
    df_results = pd.DataFrame(dic_results)

    strats = old_strats # Retomando o global original.
    columns_order = ['beta']+strats
    df_results = df_results[columns_order]

    return df_results,tab1


def get_df_resample(df_orig):
    sample_id = df_orig['id']
    N = df_orig['id'].count()
    # Sorteia N números do range(N) (i.e. do conj. {0,...,N}).
    array_elements = np.random.choice(N, N)
    ''' Extraindo da amostra inicial os termos referentes aos N í­ndices 
    sorteados acima; i.e., criando uma subamostra, também de tamanho N.'''
    resample_id = sample_id[array_elements]
    dic_resample_id = {'id':resample_id}
    df_resample_id = pd.DataFrame(dic_resample_id)
    df_resample = df_resample_id.merge(df_orig, on=['id'], how='inner')
    return df_resample


def bootstrap_MLE(M,cd_treatments):
    global df

    dic_results = {}
    for treat in cd_treatments:
        print('****TREATMENT: %s' %(treat))
        filt_treat = (df_alltreats['treatment']==treat)
        df_orig = df_alltreats[filt_treat].reset_index(drop=True)
        for m in range(M):
            if (M+1)%10==0: print(M+1)
            df = get_df_resample(df_orig)
            dic_results_row = MLE(strats,beta_guess,cte_phi_guess,bnd_beta,bnds_phi)
            dic_results_row['treatment'] = [treat]
            dic_results_row['m'] = [m+1]
            dic_results = merge_dics(dic_results,dic_results_row)
    results = pd.DataFrame(dic_results)
    return results


def get_covar(cd_treatments):
    dic_cov_mtxs = {}
    dic_std_coef = {}
    mtx_names = ['cov_mtx_treat'+str(treat) for treat in cd_treatments]
    std_coef_names = ['std_coef_treat'+str(treat) for treat in cd_treatments]
    for treat in cd_treatments:
        filt = (bootstrap_results['treatment']==treat)

        mtx_names[treat-1] = bootstrap_results.loc[filt,main_strats].cov()
        dic_cov_mtxs[treat] = mtx_names[treat-1]

        std_coef_names[treat-1] = bootstrap_results.loc[filt,main_strats].std()
        dic_std_coef[treat] = std_coef_names[treat-1]
    df_std_coef = pd.DataFrame(dic_std_coef).transpose()
    return dic_cov_mtxs, df_std_coef 


def get_Wald_pv(cd_treatments):
    dic_Wald_p_value = {strat:[] for strat in main_strats}
    dic_Wald_p_value['treatment'] = []
    for treat in cd_treatments:
        dic_Wald_p_value['treatment'] += [treat]
        prevalences = tab_prev.loc[treat-1,main_strats]
        coefs = tab_coef.loc[treat-1,main_strats]
        var_coefs = tab_std_coef.loc[treat,main_strats]**2
        for i in range(len(prevalences)):
            W2_strat = (coefs[i]-prevalences[i])**2/var_coefs[i]
            p_value = 1 - chi2.cdf(W2_strat, 1)
            dic_Wald_p_value[main_strats[i]] += [p_value]
    df_Wald_p_value = pd.DataFrame(dic_Wald_p_value)
    return df_Wald_p_value


if __name__=="__main__":

    ''' Preencher este bloco. Com exceção do arquivo da base, todos os parâmetros são 
    retirados daqui.'''
    strats = ['ac','tft','wsls','grim','stft','ad']
    # main_strats p/ calc. das matrizes de covariância (coef. de wsls é sempre quase zero).
    main_strats = ['ac','tft','grim','stft','ad'] 
    cd_treatments = range(1,6) # correspondem aos códigos dos tratamentos na base.
    beta_guess = 0.95
    cte_phi_guess = 1/len(strats) # i.e., proporções iguais.
    bnd_beta = (0.5,1)
    bnds_phi = (0,1)
    M = 1000 # Número de iterações para o bootstrap.


    '''[Ref. à Tab5 - estimação]
      Usa a base comentada no tratamento prévio, que inclui todas as estratégias 
    e retira o primeiro quarto dos dados. Tenta estimar taxas de cooperação a partir  
    das seis estratégias principais conforme o modelo apresentado. Os coeficientes são 
    os phi de cada estratégia e o beta.'''

    df = df_alltreats
    ''' O primeiro return de MLE_treats é uma reestimação excluindo os coeficientes 
    zero (ou praticamente zero) da 1º estimação. Isso serve à Tab4. Estamos interessados 
    agora no segundo return, que é simplemente a estimação com as seis estratégias.'''
    reestimate,tab_coef = MLE_treats(cd_treatments)
    tab_coef.to_excel('Tab5_coefficients.xlsx')


    '''[Ref. à Tab5 - testes]
      A matrix de variância-covariância dos coeficientes estimados é calculada por 
    um bootstrap de 1000 re-samples. A matriz então é usada para testar a igualdade 
    estatística dos coeficientes estimados com a prevalência observada das estratégias.'''

    # Demora bastante para rodar. Para testar, tirar os '#' e escolher um M pequeno.
    # bootstrap_results = bootstrap_MLE(M,cd_treatments)
    # bootstrap_results.to_csv('bootstrap_results.csv',sep=';',index=False)

    bootstrap_results = pd.read_csv('bootstrap_results.csv',sep=';')
    main_strats = ['ac','tft','grim','stft','ad'] # Sem wsls, cujo coef. é sempre zero.

    cov_mtxs,tab_std_coef = get_covar(cd_treatments)
    tab_std_coef.to_excel('Tab5_desvio_padrao_coef.xlsx')
    print(tab_std_coef)


    # Testes de Wald univariados, uma para cada estratégia, para todos os tratamentos.

    tab_prev = pd.read_excel('Tab5_prevalence.xlsx')  # Fora gerado pela rotina acima.

    df_Wald_p_value = get_Wald_pv(cd_treatments)
    df_Wald_p_value.to_excel('Wald_p_values.xlsx')
...