BitterBuster: A virtual world for understanding naturalistic patch foraging.
The 1 match
- [1] § Results › Experiment 1: Solo foraging ↔ Scripts/AnalysisScripts/AnalysisControl.ipynb, lines 1833 to 1957 · score 0.58 · confidence interval, marginal reward rate, patch visits, histogram, exit · Report a mismatch
Paper
Loaded from Europe PMC by your browser; OSCR does not store it: doi.org · Europe PMC
The paper is loaded when this pane is shown.
The authors' code
Jupyter notebook · 2,002 lines · 71 KB · no license · 1 match
- # %%
- import numpy as np
- import sys
- from scipy.stats import norm
- sys.path.append('../ProcessingScripts/')
- #from logs_to_dfs import logs_to_dfs
- import os
- import pandas as pd
- import matplotlib.pyplot as plt
- from matplotlib.collections import LineCollection, PolyCollection
- import matplotlib.colors as mcolors
- from matplotlib.lines import Line2D
- import seaborn as sns
- import math
- from statsmodels.stats.anova import AnovaRM
- import statsmodels.formula.api as smf
- from scipy.stats import spearmanr
- from scipy.stats import t as t_dist
- import matplotlib.gridspec as gridspec
- cwd = os.getcwd()
- tier_path = os.path.dirname(os.path.dirname(cwd))
- analysis_dir = os.path.join(tier_path,'Data','analysisData')
- getting_codes = os.listdir(analysis_dir)
- subject_pair_codes = [f for f in getting_codes if not f.startswith('.')] ## no hidden dirs
- # %% [markdown]
- # # Functions
- # %%
- nameSpaces = {
- 'PatchEfficiency': 'nb',
- 'EffectiveReward': '$R_e$',
- 'TimeOnPatchAvg': 'Time-On-Patch Average',
- 'VisitedNeighborhoodDepletion': 'First Visit Depletion',
- 'DPrime': 'd\'',
- 'HarvestEfficiency': 'hb'
- }
- # %%
- def makeRoundAverageGraphV2(allData, measureName, ylimit = None):
- rounds = [0, 1, 2]
- averages = []
- for r in allData['Round'].unique():
- tempDF = allData.loc[allData['Round'] == r]
- averages.append(math.fsum(tempDF[measureName]) / len(tempDF[measureName]))
- for id in allData['PairCode'].unique():
- temp = allData.loc[allData['PairCode'] == id]
- yvals = []
- rvals = []
- for r in range(len(temp)):
- yvals.append(temp.iloc[r][measureName])
- rvals.append(int(temp.iloc[r]['Round']) - 1)
- # print(type(temp['Round']))
- # print(type(temp[measureName]))
- sns.lineplot(
- x = rvals,
- y = yvals,
- color='k',
- linewidth=.2
- )
- print(rounds)
- print(averages)
- sns.barplot(
- x = rounds,
- y = averages,
- palette='Blues',
- alpha = .75
- )
- sns.stripplot(
- data= allData,
- x= 'Round',
- y= measureName,
- palette='Blues',
- linewidth=.4
- )
- if ylimit:
- plt.ylim(ylimit)
- plt.xlabel('Round', fontsize=22)
- if measureName in nameSpaces:
- plt.ylabel(nameSpaces[measureName], fontsize=22)
- else:
- plt.ylabel(measureName, fontsize=22)
- plt.xticks(rounds, ['1 (20min)', '2 (10min)', '3 (10min)'], fontsize=14)
- plt.yticks(fontsize=14)
- sns.despine()
- #plt.show()
- plt.savefig(os.path.join(tier_path, 'Output', f'F2_{measureName}_RoundDistributions.png'))
- print(averages)
- # %%
- def makeMeasureAvgsAgainstLockout(allData, measureName, ylimit = None):
- rounds = [0, 1, 2]
- averages = {}
- for r in allData['LockoutDuration'].unique():
- tempDF = allData.loc[allData['LockoutDuration'] == r]
- averages[r] = math.fsum(tempDF[measureName]) / len(tempDF[measureName])
- # for id in allData['PairCode'].unique():
- # temp = allData.loc[allData['PairCode'] == id]
- # yvals = []
- # rvals = []
- # for r in range(len(temp)):
- # yvals.append(temp.iloc[r][measureName])
- # rvals.append(int(temp.iloc[r]['Round']) - 1)
- # # print(type(temp['Round']))
- # # print(type(temp[measureName]))
- # sns.lineplot(
- # x = rvals,
- # y = yvals,
- # color='k',
- # linewidth=.2
- # )
- allData['Round'] = allData['Round'].astype(str)
- allData['LockoutDuration'] = allData['LockoutDuration'].astype(str)
- print(rounds)
- print(averages)
- sns.barplot(
- data=allData,
- errorbar='se',
- x = 'LockoutDuration',
- y = measureName,
- order=['0', '5', '15'],
- palette='Blues',
- alpha = .75
- )
- # sns.stripplot(
- # data= allData,
- # x= 'LockoutDuration',
- # y= measureName,
- # order=['0', '5', '15'],
- # palette='Blues',
- # linewidth=.4
- # )
- if ylimit:
- plt.ylim(ylimit)
- plt.xlabel('Lockout Duration', fontsize=22)
- plt.ylabel(measureName, fontsize=16)
- #plt.xticks(['0', '5', '15'], fontsize=14)
- plt.yticks(fontsize=14)
- sns.despine()
- #plt.show()
- plt.savefig(os.path.join(tier_path, 'Output', f'F2_{measureName}_RoundDistributions.png'))
- print(averages)
- # %%
- def makeRoundAverageGraphV3(allData, measureName, ylimit = None):
- rounds = [0, 1, 2]
- averages = []
- for r in allData['Round'].unique():
- tempDF = allData.loc[allData['Round'] == r]
- averages.append(math.fsum(tempDF[measureName]) / len(tempDF[measureName]))
- # for id in allData['PairCode'].unique():
- # temp = allData.loc[allData['PairCode'] == id]
- # yvals = []
- # rvals = []
- # for r in range(len(temp)):
- # yvals.append(temp.iloc[r][measureName])
- # rvals.append(temp.iloc[r]['Round'] - 1)
- # # print(type(temp['Round']))
- # # print(type(temp[measureName]))
- # sns.lineplot(
- # x = rvals,
- # y = yvals,
- # color='k',
- # linewidth=.2
- # )
- print(rounds)
- print(averages)
- allData['Round'] = allData['Round'].astype(str)
- allData['LockoutDuration'] = allData['LockoutDuration'].astype(str)
- sns.barplot(
- data=allData,
- errorbar='se',
- x = 'Round',
- y = measureName,
- hue='LockoutDuration',
- hue_order=['0', '5', '15'],
- palette='Blues',
- alpha = .75
- )
- sns.stripplot(
- data= allData,
- x= 'Round',
- y= measureName,
- hue='LockoutDuration',
- hue_order=['0', '5', '15'],
- dodge=True,
- palette='Blues',
- linewidth=.4
- )
- if ylimit:
- plt.ylim(ylimit)
- plt.xlabel('Round', fontsize=22)
- plt.ylabel(measureName, fontsize=22)
- plt.xticks(rounds, ['1 (20min)', '2 (10min)', '3 (10min)'], fontsize=14)
- plt.yticks(fontsize=14)
- sns.despine()
- #plt.show()
- plt.legend(title='Lockout Duration', loc=(-.5, .6))
- plt.savefig(os.path.join(tier_path, 'Output', f'F2_{measureName}_RoundDistributions.png'))
- print(averages)
- # %%
- def makeRoundAverageGraphV4(allData, measureName, ylimit = None):
- '''
- Just lockout subjects
- '''
- rounds = [0, 1, 2]
- averages = []
- for r in allData['Round'].unique():
- tempDF = allData.loc[allData['Round'] == r]
- averages.append(math.fsum(tempDF[measureName]) / len(tempDF[measureName]))
- print(rounds)
- print(averages)
- allData['Round'] = allData['Round'].astype(str)
- #allData['LockoutDuration'] = allData['LockoutDuration'].astype(str)
- new = allData[allData['LockoutDuration'] != '0']
- new['LockoutDuration'] = new['LockoutDuration'].astype('category').cat.remove_unused_categories()
- sns.barplot(
- data=new,
- errorbar='se',
- x = 'Round',
- y = measureName,
- hue='LockoutDuration',
- hue_order=['5', '15'],
- palette='Blues',
- alpha = .75
- )
- if ylimit:
- plt.ylim(ylimit)
- plt.xlabel('Round', fontsize=22)
- if measureName in nameSpaces:
- plt.ylabel(nameSpaces[measureName], fontsize=16)
- else:
- plt.ylabel(measureName, fontsize=16)
- plt.xticks(rounds, ['1 (20min)', '2 (10min)', '3 (10min)'], fontsize=14)
- plt.yticks(fontsize=14)
- sns.despine()
- #plt.show()
- plt.legend(title='Lockout Duration')#, loc=(-.5, .6))
- plt.savefig(os.path.join(tier_path, 'Output', f'F2_{measureName}_RoundDistributions.png'))
- print(averages)
- # %% [markdown]
- # def makeRankAndExpectedOrderList(houseEventList, x=100):
- # output = {
- # 'visitOrder': [],
- # 'rankOrder': []
- # }
- # nextExpected = {
- # 'blue': 1, 'green': 1, 'pink': 1, 'brown': 1
- # }
- # currEntry = houseEventList[0]
- # currN = houseEventList[0]['Neighborhood']
- # count=0
- # for p in range(len(houseEventList)):
- # if houseEventList[p]['Neighborhood'] != currN:
- # currN = houseEventList[p]['Neighborhood']
- # count=0
- # output['rankOrder'] = houseEventList[p]['Rank']
- # output['visitOrder'] = nextExpected[currN]
- # nextExpected[currN] += 1
- # elif count >= x:
- # continue
- # else:
- # output['rankOrder'] = houseEventList[p]['Rank']
- # output['visitOrder'] = nextExpected[currN]
- # nextExpected[currN] += 1
- # count += 1
- # return output
- # %%
- # TODO: ADD call to add new data to this dataframe
- PerformanceDF = pd.read_csv(os.path.join(tier_path, 'Output', 'Table4Performance.csv'))
- #DemographicsDF = pd.read_csv(os.path.join(tier_path, 'Output', 'Table5Demographics.csv'))
- VisitsDF = pd.read_csv(os.path.join(tier_path, 'Output', 'Table3AllVisitEvents.csv'))
- RealExitsDF = pd.read_csv(os.path.join(tier_path, 'Output', 'Table5RealExits.csv'))
- # %%
- test = VisitsDF[VisitsDF['Round'] == 1]
- aaDepletesRound1 = test['ParticipantID'].value_counts()
- aaDepletesRound1 = aaDepletesRound1 % 43
- test = VisitsDF[VisitsDF['Round'] == 2]
- aaDepletesRound2 = test['ParticipantID'].value_counts()
- aaDepletesRound2 = aaDepletesRound2 % 43
- test = VisitsDF[VisitsDF['Round'] == 3]
- aaDepletesRound3 = test['ParticipantID'].value_counts()
- aaDepletesRound3 = aaDepletesRound3 % 43
- # %% [markdown]
- # # Outputs
- # %%
- makeRoundAverageGraphV2(PerformanceDF, 'PatchEfficiency', [0,1])
- # %% [markdown]
- # makeMeasureAvgsAgainstLockout(PerformanceDF, 'VisitedNeighborhoodDepletion', [.8, 1.0])
- # %%
- makeRoundAverageGraphV2(PerformanceDF, 'EffectiveReward')
- # %% [markdown]
- # makeRoundAverageGraphV4(PerformanceDF, 'VisitedNeighborhoodDepletion', [.8, 1.0])
- # %%
- makeRoundAverageGraphV2(PerformanceDF, 'DPrime')
- # %%
- makeRoundAverageGraphV2(PerformanceDF, 'HarvestEfficiency')
- # %%
- makeRoundAverageGraphV3(PerformanceDF, 'HitsPerMinute')
- # %%
- makeMeasureAvgsAgainstLockout(PerformanceDF, 'TimeOnPatchAvg')
- # %%
- makeRoundAverageGraphV3(PerformanceDF, 'EffectiveReward')
- # %%
- makeRoundAverageGraphV3(PerformanceDF, 'RewardRate')
- # %% [markdown]
- # ## Stats
- # %%
- outMixedLM = smf.mixedlm('EffectiveReward ~ Round + DPrime * PatchEfficiency', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('HitsPerMinute ~ Round + DPrime * PatchEfficiency', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x)
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('HitsPerMinute ~ Round + EffectiveRewardBeforeAbandon * DPrime + PatchEfficiency', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x)
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('EffectiveReward ~ Round + HarvestEfficiency', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('EffectiveReward ~ PatchEfficiency', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('RewardRate ~ Round + DPrime * PatchEfficiency', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- #print(FittedOut.summary())
- # %%
- outMixedLM = smf.mixedlm('EffectiveReward ~ Round', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('HarvestEfficiency ~ Round', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('PatchEfficiency ~ Round', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x.as_latex())
- # %%
- outMixedLM = smf.mixedlm('DPrime ~ Round', data= PerformanceDF, groups=PerformanceDF['PairCode'])
- FittedOut = outMixedLM.fit()
- x = FittedOut.summary()
- #print(x.as_latex())
- # %% [markdown]
- # # To Be Redone
- # %%
- GraphingToHitsPerMin, subplots = plt.subplots(2,3, figsize=(15,10))
- rounds = ['Round1', 'Round2', 'Round3']
- GraphingToHitsPerMin.suptitle('')
- for s in range(len(rounds)):
- print(s)
- tempForSpear = PerformanceDF.loc[PerformanceDF['Round'] == str(s+1)]
- dSpear, dP = spearmanr(tempForSpear['DPrime'], tempForSpear['HitsPerMinute'])
- print(f'ROUND{s+1} Spearman R for hits-per-minute and d\' is:\n {dSpear}, p-value: {dP}\n\n')
- peSpear, peP = spearmanr(tempForSpear['PatchEfficiency'], tempForSpear['HitsPerMinute'])
- print(f'ROUND{s+1} Spearman R for hits-per-minute and PE is:\n {peSpear}, p-value: {peP}\n\n')
- sns.scatterplot(
- ax= subplots[0, s],
- data= tempForSpear,
- y= f'PatchEfficiency',
- x= f'HitsPerMinute',
- hue= f'DPrime',
- hue_norm=(0, 5),
- palette='Purples',
- linewidth=.5,
- edgecolor='k'
- )
- sns.scatterplot(
- ax= subplots[1, s],
- data= tempForSpear,
- y= f'DPrime',
- x= f'HitsPerMinute',
- hue= f'PatchEfficiency',
- hue_norm=(0.0, 1.0),
- palette='Purples',
- linewidth=.5,
- edgecolor='k'
- )
- subplots[0,s].set_xlim([0, 18])
- subplots[0,s].set_ylim([0, 1])
- subplots[0,s].set_xlabel('')
- subplots[0,s].set_ylabel('', fontsize=12)
- subplots[0,s].set_title(f'Round {s+1}', fontsize=22)
- #subplots[0,s].text(0, .85, s=f' r: {peSpear}\n p: {peP}')
- subplots[1,s].set_xlim([0, 18])
- subplots[1,s].set_ylim([-1, 5])
- subplots[1,s].set_xlabel('Hits Per Minute', fontsize=16)
- subplots[1,s].set_ylabel('', fontsize=12)
- #subplots[1,s].set_title(f'\nRound {s+1}', fontsize=16)
- #subplots[1,s].text(0, 5, s=f' r: {dSpear}\n p: {dP}')
- if s != 1:
- subplots[0,s].get_legend().set_visible(False)
- subplots[1,s].get_legend().set_visible(False)
- else:
- subplots[0,s].legend(title='d\'', loc=(-.5, .6))
- subplots[1,s].legend(title='nb', loc=(-.5, .6))
- subplots[0,0].set_ylabel('Neighborhood Bias', fontsize=16)
- subplots[1,0].set_ylabel('d\'', fontsize=16)
- sns.despine()
- #fontsize 18
- # %%
- GraphingToHitsPerMin, subplots = plt.subplots(2,3, figsize=(15,10))
- rounds = ['Round1', 'Round2', 'Round3']
- GraphingToHitsPerMin.suptitle('')
- for s in range(len(rounds)):
- tempForSpear = PerformanceDF.loc[PerformanceDF['Round'] == str(s+1)]
- dSpear, dP = spearmanr(tempForSpear['DPrime'], tempForSpear['EffectiveReward'])
- print(f'ROUND{s+1} Spearman R for effective reward and d\' is:\n {dSpear}, p-value: {dP}\n\n')
- peSpear, peP = spearmanr(tempForSpear['PatchEfficiency'], tempForSpear['EffectiveReward'])
- print(f'ROUND{s+1} Spearman R for effective reward and PE is:\n {peSpear}, p-value: {peP}\n\n')
- sns.scatterplot(
- ax= subplots[0, s],
- data= tempForSpear,
- x= f'PatchEfficiency',
- y= f'EffectiveReward',
- hue= f'DPrime',
- hue_norm=(0, 5),
- palette='Purples',
- linewidth=.5,
- edgecolor='k'
- )
- sns.scatterplot(
- ax= subplots[1, s],
- data= tempForSpear,
- x= f'DPrime',
- y= f'EffectiveReward',
- hue= f'PatchEfficiency',
- hue_norm=(0.0, 1.0),
- palette='Purples',
- linewidth=.5,
- edgecolor='k'
- )
- subplots[0,s].set_ylim([0, 0.25])
- subplots[0,s].set_xlim([0, 1])
- subplots[0,s].set_ylabel('')
- subplots[0,s].set_xlabel('nb', fontsize=16)
- subplots[0,s].set_title(f'Round {s+1}', fontsize=22)
- #subplots[0,s].text(0, .85, s=f' r: {peSpear}\n p: {peP}')
- subplots[1,s].set_ylim([0, 0.25])
- subplots[1,s].set_xlim([-1, 5])
- subplots[1,s].set_xlabel('d\'', fontsize=16)
- subplots[1,s].set_ylabel('', fontsize=12)
- #subplots[1,s].set_title(f'\nRound {s+1}', fontsize=16)
- #subplots[1,s].text(0, 5, s=f' r: {dSpear}\n p: {dP}')
- if s != 1:
- subplots[0,s].get_legend().set_visible(False)
- subplots[1,s].get_legend().set_visible(False)
- else:
- subplots[0,s].legend(title='d\'', loc=(-.5, .6))
- subplots[1,s].legend(title='nb', loc=(-.5, .6))
- subplots[0,0].set_ylabel('$R_e$', fontsize=16)
- subplots[1,0].set_ylabel('$R_e$', fontsize=16)
- sns.despine()
- #fontsize 18
- plt.savefig(os.path.join(tier_path, 'Output', 'F3_EffectiveReward_by_Performance.png'))
- # %%
- stepSize = 2.0
- endTimes = {1: 1200, 2: 600, 3: 600}
- # sort chronologically
- VisitsDF = VisitsDF.sort_values(['Round', 'ParticipantID', 'TimeIn'])
- processed_dfs = []
- for roundN in [1, 2, 3]:
- visitsForRound = VisitsDF[VisitsDF['Round'] == roundN]
- for pID in visitsForRound['ParticipantID'].unique():
- df_p = visitsForRound[visitsForRound['ParticipantID'] == pID].copy()
- if df_p.empty:
- continue
- # Build Timers & Neighborhood Visit Segments
- df_p['NbhChange'] = (df_p['NeighborhoodID'] != df_p['NeighborhoodID'].shift(1))
- df_p['NbhSegmentID'] = df_p['NbhChange'].cumsum().astype(int)
- # if df_p['TimeIn'].iloc[0] != np.nan:
- # df_p['firstIn'] = 0
- # else:
- # df_p['firstIn'] = np.nan
- df_p['firstIn'] = df_p['TimeIn'].iloc[0]
- df_p['nbhIn'] = df_p.groupby('NbhSegmentID')['TimeIn'].transform('first')
- # Build Cumulative Rewards based on TimeOut
- df_p_out = df_p.sort_values('TimeOut').copy()
- df_p_out['CumTotalReward'] = df_p_out['RewardCollected'].cumsum()
- df_p_out['CumNbhReward'] = df_p_out.groupby('NbhSegmentID')['RewardCollected'].cumsum()
- # 2 second interval
- grid_df = pd.DataFrame({
- 'RegularTimestamp': np.arange(0, endTimes[roundN] + stepSize, stepSize)
- })
- # 1. Merge active timers & active NbhSegmentID based on TimeIn <= t
- resampled = pd.merge_asof(
- grid_df,
- df_p[['TimeIn', 'firstIn', 'nbhIn', 'NbhSegmentID']],
- left_on='RegularTimestamp',
- right_on='TimeIn',
- direction='backward'
- )
- # 2. Merge Total Cumulative Reward achieved by TimeOut <= t
- resampled = pd.merge_asof(
- resampled,
- df_p_out[['TimeOut', 'CumTotalReward']],
- left_on='RegularTimestamp',
- right_on='TimeOut',
- direction='backward'
- )
- # 3. Merge Current Neighborhood Cumulative Reward achieved by TimeOut <= t
- resampled['NbhSegmentID_tmp'] = resampled['NbhSegmentID'].fillna(-1).astype(int)
- df_p_out['NbhSegmentID_tmp'] = df_p_out['NbhSegmentID'].astype(int)
- resampled = pd.merge_asof(
- resampled,
- df_p_out[['TimeOut', 'NbhSegmentID_tmp', 'CumNbhReward']],
- left_on='RegularTimestamp',
- right_on='TimeOut',
- by='NbhSegmentID_tmp',
- direction='backward'
- )
- resampled.drop(columns=['NbhSegmentID_tmp'], inplace=True)
- resampled['CumTotalReward'] = resampled['CumTotalReward'].fillna(0)
- resampled['CumNbhReward'] = resampled['CumNbhReward'].fillna(0)
- t = resampled['RegularTimestamp']
- dt_total = t# - resampled['firstIn']
- dt_nbh = t - resampled['nbhIn']
- # Calculate rates only when elapsed time > 0
- resampled['RewardRateTotal'] = np.where(dt_total > 0, resampled['CumTotalReward'] / dt_total, np.nan)
- resampled['RewardRateCurrent'] = np.where(dt_nbh > 0, resampled['CumNbhReward'] / dt_nbh, np.nan)
- resampled['VisitNumber'] = resampled['NbhSegmentID'] # visit tracking columns for multi-panel plotting
- resampled['TimeInNbh'] = np.where(dt_nbh >= 0, dt_nbh, np.nan)
- resampled['ID'] = pID
- resampled['Round'] = roundN
- processed_dfs.append(resampled[[
- 'ID', 'Round', 'RegularTimestamp',
- 'RewardRateTotal', 'RewardRateCurrent',
- 'VisitNumber', 'TimeInNbh', 'nbhIn', 'firstIn'
- ]])
- resampledDF = pd.concat(processed_dfs, ignore_index=True)
- # %%
- '''
- Bypass this if not using cumulative expected reward rate
- '''
- stepSize = 2.0
- endTimes = {1: 1200, 2: 600, 3: 600}
- # Global offsets to map local round time to continuous experiment time
- round_offsets = {1: 0, 2: 1200, 3: 1800}
- # Sort chronologically by Participant and Round
- VisitsDF = VisitsDF.sort_values(['ParticipantID', 'Round', 'TimeIn']).copy()
- # Add global timestamps for experiment-wide matching
- VisitsDF['GlobalTimeIn'] = VisitsDF['TimeIn'] + VisitsDF['Round'].map(round_offsets)
- VisitsDF['GlobalTimeOut'] = VisitsDF['TimeOut'] + VisitsDF['Round'].map(round_offsets)
- processed_dfs = []
- # Process per participant across all rounds
- for pID in VisitsDF['ParticipantID'].unique():
- df_p = VisitsDF[VisitsDF['ParticipantID'] == pID].copy()
- if df_p.empty:
- continue
- # 1. Global cumulative reward tracking (across all 3 rounds)
- df_p_out = df_p.sort_values('GlobalTimeOut').copy()
- df_p_out['CumTotalReward'] = df_p_out['RewardCollected'].cumsum()
- # 2. Process each round with reset round-level timers
- for roundN in [1, 2, 3]:
- df_p_round = df_p[df_p['Round'] == roundN].copy()
- if df_p_round.empty:
- continue
- # Local Neighborhood Tracking (resets/tracks for current round)
- df_p_round['NbhChange'] = (df_p_round['NeighborhoodID'] != df_p_round['NeighborhoodID'].shift(1))
- df_p_round['NbhSegmentID'] = df_p_round['NbhChange'].cumsum().astype(int)
- #df_p_round['firstIn'] = 0 if pd.notna(df_p_round['TimeIn'].iloc[0]) else np.nan
- df_p_round['firstIn'] = df_p_round['TimeIn'].iloc[0]
- df_p_round['nbhIn'] = df_p_round.groupby('NbhSegmentID')['TimeIn'].transform('first')
- # Local Cumulative Neighborhood Reward
- df_p_round_out = df_p_round.sort_values('TimeOut').copy()
- df_p_round_out['CumNbhReward'] = df_p_round_out.groupby('NbhSegmentID')['RewardCollected'].cumsum()
- # Regular local grid starting at 0 for plotting
- grid_df = pd.DataFrame({
- 'RegularTimestamp': np.arange(0, endTimes[roundN] + stepSize, stepSize)
- })
- # Global lookup timestamp for cross-round reward lookup
- grid_df['GlobalTimestamp'] = grid_df['RegularTimestamp'] + round_offsets[roundN]
- # Merge active local timers & NbhSegmentID
- resampled = pd.merge_asof(
- grid_df,
- df_p_round[['TimeIn', 'firstIn', 'nbhIn', 'NbhSegmentID']],
- left_on='RegularTimestamp',
- right_on='TimeIn',
- direction='backward'
- )
- # Merge ALL-TIME Cumulative Reward up to global time point
- resampled = pd.merge_asof(
- resampled,
- df_p_out[['GlobalTimeOut', 'CumTotalReward']],
- left_on='GlobalTimestamp',
- right_on='GlobalTimeOut',
- direction='backward'
- )
- # Merge CURRENT Neighborhood Cumulative Reward up to local time point
- resampled['NbhSegmentID_tmp'] = resampled['NbhSegmentID'].fillna(-1).astype(int)
- df_p_round_out['NbhSegmentID_tmp'] = df_p_round_out['NbhSegmentID'].astype(int)
- resampled = pd.merge_asof(
- resampled,
- df_p_round_out[['TimeOut', 'NbhSegmentID_tmp', 'CumNbhReward']],
- left_on='RegularTimestamp',
- right_on='TimeOut',
- by='NbhSegmentID_tmp',
- direction='backward'
- )
- resampled.drop(columns=['NbhSegmentID_tmp'], inplace=True)
- resampled['CumTotalReward'] = resampled['CumTotalReward'].fillna(0)
- resampled['CumNbhReward'] = resampled['CumNbhReward'].fillna(0)
- # Local elapsed time for round and neighborhood
- t_local = resampled['RegularTimestamp']
- dt_total_local = t_local - resampled['firstIn']
- dt_nbh_local = t_local - resampled['nbhIn']
- # Global elapsed time for cross-round total reward rate
- t_global = resampled['GlobalTimestamp']
- # Rates:
- # Total rate uses cumulative rewards over global elapsed time
- # Current rate uses local neighborhood rewards over local neighborhood time
- resampled['RewardRateTotal'] = np.where(t_global > 0, resampled['CumTotalReward'] / t_global, np.nan)
- resampled['RewardRateCurrent'] = np.where(dt_nbh_local > 0, resampled['CumNbhReward'] / dt_nbh_local, np.nan)
- resampled['VisitNumber'] = resampled['NbhSegmentID']
- resampled['TimeInNbh'] = np.where(dt_nbh_local >= 0, dt_nbh_local, np.nan)
- resampled['ID'] = pID
- resampled['Round'] = roundN
- processed_dfs.append(resampled[[
- 'ID', 'Round', 'RegularTimestamp',
- 'RewardRateTotal', 'RewardRateCurrent',
- 'VisitNumber', 'TimeInNbh', 'nbhIn', 'firstIn'
- ]])
- resampledDF = pd.concat(processed_dfs, ignore_index=True)
- # %%
- sns.set_theme(style="ticks")
- palette = {"tRR": "#1f77b4", "cRR": "#ff7f0e"}
- fig, axes = plt.subplots(3, 1, figsize=(10, 11), sharex=True, sharey=True)
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_data = resampledDF[resampledDF['Round'] == roundN]
- '''
- Plot individual participant traces
- '''
- for pID, p_data in round_data.groupby('ID'):
- ax.plot(
- p_data['RegularTimestamp'], p_data['RewardRateTotal'],
- color=palette['tRR'], alpha=0.12, linewidth=0.8, zorder=1
- )
- ax.plot(
- p_data['RegularTimestamp'], p_data['RewardRateCurrent'],
- color=palette['cRR'], alpha=0.12, linewidth=0.8, zorder=1
- )
- '''
- Calculate and plot group mean
- '''
- mean_df = round_data.groupby('RegularTimestamp')[['RewardRateTotal', 'RewardRateCurrent']].mean().reset_index()
- ax.plot(
- mean_df['RegularTimestamp'], mean_df['RewardRateTotal'],
- color=palette['tRR'], linewidth=2.5, label='Total Reward Rate (tRR)', zorder=2
- )
- ax.plot(
- mean_df['RegularTimestamp'], mean_df['RewardRateCurrent'],
- color=palette['cRR'], linewidth=2.5, label='Current Nbh Reward Rate (cRR)', zorder=2
- )
- ax.set_title(f'Round {roundN}', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Reward Rate (pts/sec)', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- if i == 0:
- ax.legend(
- loc='upper right', frameon=True, facecolor='white',
- edgecolor='none', framealpha=0.9, fontsize=10
- )
- axes[-1].set_xlabel('Time (seconds)', fontsize=11, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- palette = {"tRR": "#1f77b4", "cRR": "#ff7f0e"}
- plot_df = resampledDF.copy()
- def mask_before_first_payload(df_p):
- valid_payloads = df_p[df_p['RewardRateTotal'].notna() & (df_p['RewardRateTotal'] > 0)]
- if not valid_payloads.empty:
- first_payload_t = valid_payloads['RegularTimestamp'].iloc[0]
- df_p.loc[df_p['RegularTimestamp'] < first_payload_t, ['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- else:
- df_p[['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- return df_p
- plot_df = plot_df.groupby(['Round', 'ID'], group_keys=False).apply(mask_before_first_payload)
- fig, axes = plt.subplots(3, 1, figsize=(10, 11), sharex=True, sharey=True)
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_data = plot_df[plot_df['Round'] == roundN]
- # individual traces starting at personal first payload
- for pID, p_data in round_data.groupby('ID'):
- p_valid = p_data.dropna(subset=['RewardRateTotal', 'RewardRateCurrent'], how='all')
- if p_valid.empty:
- continue
- ax.plot(
- p_valid['RegularTimestamp'], p_valid['RewardRateTotal'],
- color=palette['tRR'], alpha=0.12, linewidth=0.8, zorder=1
- )
- ax.plot(
- p_valid['RegularTimestamp'], p_valid['RewardRateCurrent'],
- color=palette['cRR'], alpha=0.12, linewidth=0.8, zorder=1
- )
- # population mean
- mean_df = round_data.groupby('RegularTimestamp')[['RewardRateTotal', 'RewardRateCurrent']].mean().reset_index()
- mean_df = mean_df.dropna(subset=['RewardRateTotal', 'RewardRateCurrent'], how='all')
- ax.plot(
- mean_df['RegularTimestamp'], mean_df['RewardRateTotal'],
- color=palette['tRR'], linewidth=2.5, label='Total Reward Rate (tRR)', zorder=2
- )
- ax.plot(
- mean_df['RegularTimestamp'], mean_df['RewardRateCurrent'],
- color=palette['cRR'], linewidth=2.5, label='Current Nbh Reward Rate (cRR)', zorder=2
- )
- ax.set_title(f'Round {roundN}', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Reward Rate (pts/sec)', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- if i == 0:
- ax.legend(
- loc='upper right', frameon=True, facecolor='white',
- edgecolor='none', framealpha=0.9, fontsize=10
- )
- axes[-1].set_xlabel('Time (seconds)', fontsize=11, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- blue_col = "#1f77b4"
- orange_col = "#ff7f0e"
- MAX_VISITS_TO_SHOW = 5
- fig = plt.figure(figsize=(16, 10))
- outer_gs = gridspec.GridSpec(3, 1, figure=fig, hspace=0.55)
- for i, roundN in enumerate([1, 2, 3]):
- # subgrid: 1 main tRR panel + MAX_VISITS_TO_SHOW small cRR panels
- inner_gs = gridspec.GridSpecFromSubplotSpec(
- 1, 1 + MAX_VISITS_TO_SHOW,
- subplot_spec=outer_gs[i],
- width_ratios=[2.5] + [1.0] * MAX_VISITS_TO_SHOW,
- wspace=0.35
- )
- round_data = resampledDF[resampledDF['Round'] == roundN]
- '''
- MAIN PANEL (Left): Total Reward Rate (tRR) over Round Time
- '''
- ax_trr = fig.add_subplot(inner_gs[0, 0])
- # Participant background lines
- for pID, p_data in round_data.groupby('ID'):
- ax_trr.plot(p_data['RegularTimestamp'], p_data['RewardRateTotal'],
- color=blue_col, alpha=0.12, linewidth=0.8, zorder=1)
- # mean line
- mean_trr = round_data.groupby('RegularTimestamp')['RewardRateTotal'].mean().reset_index()
- ax_trr.plot(mean_trr['RegularTimestamp'], mean_trr['RewardRateTotal'],
- color=blue_col, linewidth=2.2, label='tRR (Total Rate)', zorder=2)
- ax_trr.set_title(f'Round {roundN}: Total Rate', fontsize=11, fontweight='semibold', loc='left')
- ax_trr.set_xlabel('Round Time (s)', fontsize=9)
- ax_trr.set_ylabel('Rate (reward/s)', fontsize=9)
- ax_trr.grid(True, linestyle='--', alpha=0.3)
- sns.despine(ax=ax_trr)
- if i == 0:
- ax_trr.legend(loc='upper right', frameon=True, facecolor='white', framealpha=0.9, fontsize=9)
- '''
- SUB-PANELS (Right): cRR for Visit 1, Visit 2, ..., Visit N
- '''
- for visit_idx in range(1, MAX_VISITS_TO_SHOW + 1):
- ax_crr = fig.add_subplot(inner_gs[0, visit_idx])
- visit_data = round_data[round_data['VisitNumber'] == visit_idx]
- valid_visit_data = visit_data.dropna(subset=['RewardRateCurrent'])
- n_subjects = valid_visit_data['ID'].nunique()
- if not valid_visit_data.empty:
- for pID, p_data in valid_visit_data.groupby('ID'):
- ax_crr.plot(p_data['TimeInNbh'], p_data['RewardRateCurrent'],
- color=orange_col, alpha=0.15, linewidth=0.8, zorder=1)
- mean_crr = valid_visit_data.groupby('TimeInNbh')['RewardRateCurrent'].mean().reset_index()
- ax_crr.plot(mean_crr['TimeInNbh'], mean_crr['RewardRateCurrent'],
- color=orange_col, linewidth=2.0, zorder=2)
- ax_crr.set_title(f'Visit #{visit_idx}\n(N={n_subjects})', fontsize=9, fontweight='medium')
- ax_crr.set_xlabel('In-Nbh Time (s)', fontsize=8)
- ax_crr.grid(True, linestyle='--', alpha=0.3)
- sns.despine(ax=ax_crr)
- if visit_idx > 1:
- ax_crr.set_yticklabels([])
- # global y-axis scaling
- global_ymax = max(ax.get_ylim()[1] for ax in fig.axes)
- target_ymax = max(0.5, global_ymax)
- for ax in fig.axes:
- ax.set_ylim(bottom=0, top=target_ymax)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- palette = {"tRR": "#1f77b4", "cRR": "#ff7f0e"}
- plot_df = resampledDF.copy()
- # 1. Mask timestamps prior to the first completed payload in EACH neighborhood visit
- def mask_before_first_nbh_payload(df_visit):
- valid_payloads = df_visit[df_visit['RewardRateCurrent'].notna() & (df_visit['RewardRateCurrent'] > 0)]
- if not valid_payloads.empty:
- first_payload_t = valid_payloads['TimeInNbh'].iloc[0]
- df_visit.loc[df_visit['TimeInNbh'] < first_payload_t, ['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- else:
- df_visit[['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- return df_visit
- plot_df = plot_df.groupby(['Round', 'ID', 'VisitNumber'], group_keys=False).apply(mask_before_first_nbh_payload)
- # 2. Bin TimeInNbh to integer seconds and truncate at t = 300s
- plot_df['TimeInNbh_bin'] = plot_df['TimeInNbh'].round(0)
- plot_df = plot_df[plot_df['TimeInNbh_bin'] <= 300]
- # Helper function to plot a line with dynamic opacity mapped to N
- def plot_fading_line(ax, x, y, counts, color, linewidth, zorder):
- x = np.asarray(x)
- y = np.asarray(y)
- counts = np.asarray(counts)
- n_max = counts.max()
- n_min = counts.min()
- # Scale alpha linearly between 0.20 (at N_min) and 1.00 (at N_max)
- if n_max > n_min:
- alphas = 0.20 + 0.80 * (counts - n_min) / (n_max - n_min)
- else:
- alphas = np.ones_like(counts, dtype=float)
- # Build 2D line segments between consecutive (x, y) points
- points = np.array([x, y]).T.reshape(-1, 1, 2)
- segments = np.concatenate([points[:-1], points[1:]], axis=1)
- # Calculate average RGBA color per segment
- base_rgb = mcolors.to_rgb(color)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- colors = [(*base_rgb, a) for a in seg_alphas]
- lc = LineCollection(segments, colors=colors, linewidths=linewidth, zorder=zorder)
- ax.add_collection(lc)
- fig, axes = plt.subplots(3, 1, figsize=(10, 11), sharex=True, sharey=True)
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_data = plot_df[plot_df['Round'] == roundN]
- # 3. Aggregate mean and sample size N per time bin (up to t=300)
- agg_df = (
- round_data.groupby('TimeInNbh_bin')
- .agg(
- trr_mean=('RewardRateTotal', 'mean'),
- trr_count=('RewardRateTotal', 'count'),
- crr_mean=('RewardRateCurrent', 'mean'),
- crr_count=('RewardRateCurrent', 'count')
- )
- .reset_index()
- .sort_values('TimeInNbh_bin')
- .dropna(subset=['trr_mean', 'crr_mean'], how='all')
- )
- time_vals = agg_df['TimeInNbh_bin'].values
- # 4. Plot cRR (Orange) with fading opacity
- plot_fading_line(
- ax, time_vals, agg_df['crr_mean'], agg_df['crr_count'],
- color=palette['cRR'], linewidth=2.5, zorder=2
- )
- # 5. Plot tRR (Blue) on top with fading opacity
- plot_fading_line(
- ax, time_vals, agg_df['trr_mean'], agg_df['trr_count'],
- color=palette['tRR'], linewidth=2.5, zorder=3
- )
- ax.set_title(f'Round {roundN}', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Reward Rate (reward/sec)', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- if i == 0:
- # Custom legend elements for metrics + opacity gradient indicators
- legend_elements = [
- Line2D([0], [0], color=palette['cRR'], lw=2.5, label='Patch Reward Rate (pRR)'),
- Line2D([0], [0], color=palette['tRR'], lw=2.5, label='Expected Reward Rate (eRR)'),
- Line2D([0], [0], color='gray', lw=2.5, alpha=1.0, label=r'100% Opacity ($N = N_{\mathrm{max}}$)'),
- Line2D([0], [0], color='gray', lw=2.5, alpha=0.2, label=r'20% Opacity ($N = N_{\mathrm{min}}$)'),
- ]
- ax.legend(
- handles=legend_elements,
- loc='upper right', frameon=True, facecolor='white',
- edgecolor='none', framealpha=0.9, fontsize=9.5
- )
- # Fix x-axis range to [0, 300] seconds
- for ax in axes:
- ax.set_xlim(left=0, right=300)
- # Axis formatting
- axes[-1].set_xlabel('Time Since Neighborhood Entry (seconds)', fontsize=11, labelpad=8)
- # Enforce uniform bottom=0 and a minimum y-axis ceiling across subplots
- global_ymax = max(ax.get_ylim()[1] for ax in fig.axes)
- target_ymax = max(0.5, global_ymax)
- for ax in fig.axes:
- ax.set_ylim(bottom=0, top=target_ymax)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- palette = {"tRR": "#1f77b4", "cRR": "#ff7f0e"}
- plot_df = resampledDF.copy()
- # 1. Mask timestamps prior to the first completed payload in EACH neighborhood visit
- def mask_before_first_nbh_payload(df_visit):
- valid_payloads = df_visit[df_visit['RewardRateCurrent'].notna() & (df_visit['RewardRateCurrent'] > 0)]
- if not valid_payloads.empty:
- first_payload_t = valid_payloads['TimeInNbh'].iloc[0]
- df_visit.loc[df_visit['TimeInNbh'] < first_payload_t, ['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- else:
- df_visit[['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- return df_visit
- plot_df = plot_df.groupby(['Round', 'ID', 'VisitNumber'], group_keys=False).apply(mask_before_first_nbh_payload)
- # 2. Bin TimeInNbh to integer seconds and truncate at t = 300s
- plot_df['TimeInNbh_bin'] = plot_df['TimeInNbh'].round(0)
- plot_df = plot_df[plot_df['TimeInNbh_bin'] <= 300]
- # Helper function to plot line segments with dynamic opacity mapped to N
- def plot_fading_line(ax, x, y, counts, color, linewidth, zorder):
- x = np.asarray(x)
- y = np.asarray(y)
- counts = np.asarray(counts)
- n_max = counts.max()
- n_min = counts.min()
- if n_max > n_min:
- alphas = 0.20 + 0.80 * (counts - n_min) / (n_max - n_min)
- else:
- alphas = np.ones_like(counts, dtype=float)
- points = np.array([x, y]).T.reshape(-1, 1, 2)
- segments = np.concatenate([points[:-1], points[1:]], axis=1)
- base_rgb = mcolors.to_rgb(color)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- colors = [(*base_rgb, a) for a in seg_alphas]
- lc = LineCollection(segments, colors=colors, linewidths=linewidth, zorder=zorder)
- ax.add_collection(lc)
- # Helper function to plot SEM shaded ribbon with dynamic opacity mapped to N
- def plot_fading_ribbon(ax, x, y_lower, y_upper, counts, color, zorder):
- x = np.asarray(x)
- y_lower = np.asarray(y_lower)
- y_upper = np.asarray(y_upper)
- counts = np.asarray(counts)
- n_max = counts.max()
- n_min = counts.min()
- # Scale ribbon fill opacity from 0.05 (at N_min) to 0.25 (at N_max)
- if n_max > n_min:
- alphas = 0.05 + 0.20 * (counts - n_min) / (n_max - n_min)
- else:
- alphas = np.full_like(counts, 0.25, dtype=float)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- base_rgb = mcolors.to_rgb(color)
- facecolors = [(*base_rgb, a) for a in seg_alphas]
- verts = []
- for i in range(len(x) - 1):
- verts.append([
- [x[i], y_lower[i]],
- [x[i+1], y_lower[i+1]],
- [x[i+1], y_upper[i+1]],
- [x[i], y_upper[i]]
- ])
- pc = PolyCollection(verts, facecolors=facecolors, edgecolors='none', zorder=zorder)
- ax.add_collection(pc)
- fig, axes = plt.subplots(3, 1, figsize=(10, 11), sharex=True, sharey=True)
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_data = plot_df[plot_df['Round'] == roundN]
- # 3. Aggregate mean, SEM, and sample size N per time bin (up to t=300)
- agg_df = (
- round_data.groupby('TimeInNbh_bin')
- .agg(
- trr_mean=('RewardRateTotal', 'mean'),
- trr_sem=('RewardRateTotal', 'sem'),
- trr_count=('RewardRateTotal', 'count'),
- crr_mean=('RewardRateCurrent', 'mean'),
- crr_sem=('RewardRateCurrent', 'sem'),
- crr_count=('RewardRateCurrent', 'count')
- )
- .reset_index()
- .sort_values('TimeInNbh_bin')
- .dropna(subset=['trr_mean', 'crr_mean'], how='all')
- )
- # Handle bins with single observation (SEM=NaN)
- agg_df['trr_sem'] = agg_df['trr_sem'].fillna(0)
- agg_df['crr_sem'] = agg_df['crr_sem'].fillna(0)
- time_vals = agg_df['TimeInNbh_bin'].values
- # 4. Plot cRR (Orange) SEM Ribbon + Mean Line (Underneath)
- plot_fading_ribbon(
- ax, time_vals,
- agg_df['crr_mean'] - agg_df['crr_sem'],
- agg_df['crr_mean'] + agg_df['crr_sem'],
- agg_df['crr_count'], color=palette['cRR'], zorder=1
- )
- plot_fading_line(
- ax, time_vals, agg_df['crr_mean'], agg_df['crr_count'],
- color=palette['cRR'], linewidth=2.5, zorder=2
- )
- # 5. Plot tRR (Blue) SEM Ribbon + Mean Line (On Top)
- plot_fading_ribbon(
- ax, time_vals,
- agg_df['trr_mean'] - agg_df['trr_sem'],
- agg_df['trr_mean'] + agg_df['trr_sem'],
- agg_df['trr_count'], color=palette['tRR'], zorder=3
- )
- plot_fading_line(
- ax, time_vals, agg_df['trr_mean'], agg_df['trr_count'],
- color=palette['tRR'], linewidth=2.5, zorder=4
- )
- ax.set_title(f'Round {roundN}', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Reward Rate (reward/sec)', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- if i == 0:
- legend_elements = [
- Line2D([0], [0], color=palette['cRR'], lw=2.5, label='Patch Reward Rate (pRR) ± SEM'),
- Line2D([0], [0], color=palette['tRR'], lw=2.5, label='Expected Reward Rate (eRR) ± SEM'),
- Line2D([0], [0], color='gray', lw=2.5, alpha=1.0, label=r'100% Opacity ($N = N_{\mathrm{max}}$)'),
- Line2D([0], [0], color='gray', lw=2.5, alpha=0.2, label=r'20% Opacity ($N = N_{\mathrm{min}}$)'),
- ]
- ax.legend(
- handles=legend_elements,
- loc='upper right', frameon=True, facecolor='white',
- edgecolor='none', framealpha=0.9, fontsize=9.5
- )
- # Fix x-axis range to [0, 300] seconds
- for ax in axes:
- ax.set_xlim(left=0, right=300)
- # Axis formatting
- axes[-1].set_xlabel('Time Since Neighborhood Entry (seconds)', fontsize=11, labelpad=8)
- # Enforce uniform bottom=0 and a minimum y-axis ceiling across subplots
- global_ymax = max(ax.get_ylim()[1] for ax in fig.axes)
- target_ymax = max(0.5, global_ymax)
- for ax in fig.axes:
- ax.set_ylim(bottom=0, top=target_ymax)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- palette = {"tRR": "#1f77b4", "cRR": "#ff7f0e"}
- plot_df = resampledDF.copy()
- # 1. Mask timestamps prior to the first completed payload in EACH neighborhood visit
- def mask_before_first_nbh_payload(df_visit):
- valid_payloads = df_visit[df_visit['RewardRateCurrent'].notna() & (df_visit['RewardRateCurrent'] > 0)]
- if not valid_payloads.empty:
- first_payload_t = valid_payloads['TimeInNbh'].iloc[0]
- df_visit.loc[df_visit['TimeInNbh'] < first_payload_t, ['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- else:
- df_visit[['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- return df_visit
- plot_df = plot_df.groupby(['Round', 'ID', 'VisitNumber'], group_keys=False).apply(mask_before_first_nbh_payload)
- # 2. Bin TimeInNbh to integer seconds and truncate at t = 300s
- plot_df['TimeInNbh_bin'] = plot_df['TimeInNbh'].round(0)
- plot_df = plot_df[plot_df['TimeInNbh_bin'] <= 300]
- # Helper function to plot a line with dynamic opacity mapped to N
- def plot_fading_line(ax, x, y, counts, color, linewidth, zorder):
- x = np.asarray(x)
- y = np.asarray(y)
- counts = np.asarray(counts)
- n_max = counts.max()
- n_min = counts.min()
- # Scale alpha linearly between 0.20 (at N_min) and 1.00 (at N_max)
- if n_max > n_min:
- alphas = 0.20 + 0.80 * (counts - n_min) / (n_max - n_min)
- else:
- alphas = np.ones_like(counts, dtype=float)
- # Build 2D line segments between consecutive (x, y) points
- points = np.array([x, y]).T.reshape(-1, 1, 2)
- segments = np.concatenate([points[:-1], points[1:]], axis=1)
- # Calculate average RGBA color per segment
- base_rgb = mcolors.to_rgb(color)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- colors = [(*base_rgb, a) for a in seg_alphas]
- lc = LineCollection(segments, colors=colors, linewidths=linewidth, zorder=zorder)
- ax.add_collection(lc)
- fig, axes = plt.subplots(3, 1, figsize=(10, 11), sharex=True, sharey=True)
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_data = plot_df[plot_df['Round'] == roundN]
- # 3. Individual visit traces underneath (linewidth=0.4, alpha=0.12, zorder=1)
- for (pID, vNum), visit_data in round_data.groupby(['ID', 'VisitNumber']):
- v_valid = visit_data.dropna(subset=['RewardRateTotal', 'RewardRateCurrent'], how='all')
- if v_valid.empty:
- continue
- # Orange (cRR) individual visit trace
- ax.plot(
- v_valid['TimeInNbh_bin'], v_valid['RewardRateCurrent'],
- color=palette['cRR'], alpha=0.12, linewidth=0.4, zorder=1
- )
- # Blue (tRR) individual visit trace
- ax.plot(
- v_valid['TimeInNbh_bin'], v_valid['RewardRateTotal'],
- color=palette['tRR'], alpha=0.12, linewidth=0.4, zorder=1
- )
- # 4. Aggregate mean and sample size N per time bin (up to t=300)
- agg_df = (
- round_data.groupby('TimeInNbh_bin')
- .agg(
- trr_mean=('RewardRateTotal', 'mean'),
- trr_count=('RewardRateTotal', 'count'),
- crr_mean=('RewardRateCurrent', 'mean'),
- crr_count=('RewardRateCurrent', 'count')
- )
- .reset_index()
- .sort_values('TimeInNbh_bin')
- .dropna(subset=['trr_mean', 'crr_mean'], how='all')
- )
- time_vals = agg_df['TimeInNbh_bin'].values
- # 5. Plot cRR (Orange) mean with fading opacity
- plot_fading_line(
- ax, time_vals, agg_df['crr_mean'], agg_df['crr_count'],
- color=palette['cRR'], linewidth=2.5, zorder=2
- )
- # 6. Plot tRR (Blue) mean on top with fading opacity
- plot_fading_line(
- ax, time_vals, agg_df['trr_mean'], agg_df['trr_count'],
- color=palette['tRR'], linewidth=2.5, zorder=3
- )
- ax.set_title(f'Round {roundN}', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Reward Rate (reward/sec)', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- if i == 0:
- # Custom legend elements
- legend_elements = [
- Line2D([0], [0], color=palette['cRR'], lw=2.5, label='Patch Reward Rate (pRR)'),
- Line2D([0], [0], color=palette['tRR'], lw=2.5, label='Expected Reward Rate (eRR)'),
- Line2D([0], [0], color='gray', lw=0.8, alpha=0.5, label='Individual Visit Traces'),
- Line2D([0], [0], color='gray', lw=2.5, alpha=1.0, label=r'100% Opacity ($N = N_{\mathrm{max}}$)'),
- Line2D([0], [0], color='gray', lw=2.5, alpha=0.2, label=r'20% Opacity ($N = N_{\mathrm{min}}$)'),
- ]
- ax.legend(
- handles=legend_elements,
- loc='upper right', frameon=True, facecolor='white',
- edgecolor='none', framealpha=0.9, fontsize=9.5
- )
- # Fix x-axis range to [0, 300] seconds
- for ax in axes:
- ax.set_xlim(left=0, right=300)
- # Axis formatting
- axes[-1].set_xlabel('Time Since Neighborhood Entry (seconds)', fontsize=11, labelpad=8)
- # Enforce uniform bottom=0 and a minimum y-axis ceiling across subplots
- global_ymax = max(ax.get_ylim()[1] for ax in fig.axes)
- target_ymax = max(0.5, global_ymax)
- for ax in fig.axes:
- ax.set_ylim(bottom=0, top=target_ymax)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- # --- 1. Blue Color Gradient across Rounds ---
- blues = sns.color_palette("Blues", n_colors=6)
- round_colors = [blues[2], blues[4], blues[5]] # Round 1: Light, Round 2: Medium, Round 3: Dark
- # --- 2. Vectorized Data Processing ---
- plot_df = resampledDF.copy()
- # Find time of first valid payload per visit group
- valid_mask = plot_df['RewardRateCurrent'] > 0
- plot_df['first_payload_t'] = plot_df.where(valid_mask).groupby(['Round', 'ID', 'VisitNumber'])['TimeInNbh'].transform('min')
- # Mask timestamps prior to first payload (or visits with no payload)
- invalid_mask = (plot_df['TimeInNbh'] < plot_df['first_payload_t']) | plot_df['first_payload_t'].isna()
- plot_df.loc[invalid_mask, ['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- # Calculate Marginal Reward Rate (Patch RR - Expected RR) and filter binned time (<= 200s)
- plot_df['RR_Diff'] = plot_df['RewardRateCurrent'] - plot_df['RewardRateTotal']
- plot_df['TimeInNbh_bin'] = plot_df['TimeInNbh'].round(0)
- plot_df = plot_df[plot_df['TimeInNbh_bin'] <= 400]
- # --- 3. Compact Visual Helpers ---
- def get_alphas(counts, min_val, max_val):
- """Rescale counts to alpha range [min_val, max_val]."""
- n_min, n_max = counts.min(), counts.max()
- if n_max > n_min:
- return min_val + (max_val - min_val) * (counts - n_min) / (n_max - n_min)
- return np.full_like(counts, max_val, dtype=float)
- def plot_fading_line(ax, x, y, counts, color, linewidth=2.5, zorder=2):
- alphas = get_alphas(counts, 0.20, 1.00)
- pts = np.array([x, y]).T.reshape(-1, 1, 2)
- segments = np.concatenate([pts[:-1], pts[1:]], axis=1)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- colors = [(*mcolors.to_rgb(color), a) for a in seg_alphas]
- ax.add_collection(LineCollection(segments, colors=colors, linewidths=linewidth, zorder=zorder))
- def plot_fading_ribbon(ax, x, y_lower, y_upper, counts, color, zorder=1):
- alphas = get_alphas(counts, 0.05, 0.25)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- base_rgb = mcolors.to_rgb(color)
- verts = [
- [[x[i], y_lower[i]], [x[i+1], y_lower[i+1]],
- [x[i+1], y_upper[i+1]], [x[i], y_upper[i]]]
- for i in range(len(x) - 1)
- ]
- ax.add_collection(PolyCollection(verts, facecolors=[(*base_rgb, a) for a in seg_alphas], edgecolors='none', zorder=zorder))
- # --- 4. Figure Plotting Loop ---
- fig, axes = plt.subplots(3, 1, figsize=(10, 11), sharex=True, sharey=True)
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- current_color = round_colors[i] # Progressive blue per round
- # Aggregate mean, sample size N, and compute 95% CI
- agg_df = (
- plot_df[plot_df['Round'] == roundN]
- .groupby('TimeInNbh_bin')['RR_Diff']
- .agg(mean='mean', sem='sem', count='count')
- .dropna(subset=['mean'])
- .reset_index()
- )
- # Calculate 95% Confidence Interval bound using scipy.stats.t
- t_crit = t_dist.ppf(0.975, df=np.maximum(agg_df['count'] - 1, 1))
- agg_df['ci95'] = (t_crit * agg_df['sem']).fillna(0)
- x = agg_df['TimeInNbh_bin'].values
- y_mean = agg_df['mean'].values
- ci = agg_df['ci95'].values
- counts = agg_df['count'].values
- # Plot elements
- ax.axhline(0, color='gray', linestyle='--', linewidth=1.2, alpha=0.7, zorder=0)
- plot_fading_ribbon(ax, x, y_mean - ci, y_mean + ci, counts, color=current_color, zorder=1)
- plot_fading_line(ax, x, y_mean, counts, color=current_color, linewidth=2.5, zorder=2)
- # Subplot Formatting
- ax.set_ylim(-0.2, 0.2)
- ax.set_xlim(0, 400)
- ax.set_title(f'Round {roundN}', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Marginal Reward Rate (reward/sec)\n[Patch RR - Expected RR]', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- axes[-1].set_xlabel('Time Since Neighborhood Entry (seconds)', fontsize=11, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- # --- 1. Calculate Marginal RR ---
- resampled_df = resampledDF.copy()
- resampled_df['Marginal_RR'] = resampled_df['RewardRateCurrent'] - resampled_df['RewardRateTotal']
- exits_df = RealExitsDF.copy().rename(columns={'ParticipantID': 'ID'})
- # --- 2. Align Key Data Types Across DataFrames ---
- merge_by = ['Round', 'ID']
- if 'VisitNumber' in exits_df.columns and 'VisitNumber' in resampled_df.columns:
- merge_by.append('VisitNumber')
- for col in merge_by:
- exits_df[col] = exits_df[col].astype(int)
- resampled_df[col] = resampled_df[col].astype(int)
- # --- 3. Sort & Merge on Time Keys ---
- resampled_sorted = resampled_df.dropna(subset=['Marginal_RR']).sort_values('RegularTimestamp')
- exits_sorted = exits_df.sort_values('TimeStamp')
- exit_data = pd.merge_asof(
- exits_sorted,
- resampled_sorted[merge_by + ['RegularTimestamp', 'Marginal_RR']],
- left_on='TimeStamp', # Time column in RealExitsDF
- right_on='RegularTimestamp', # Time column in resampledDF
- by=merge_by,
- direction='nearest' # Matches to nearest floating timestamp
- )
- # Diagnostic verification
- print(f"Total Exit Events in RealExitsDF: {len(RealExitsDF)}")
- print(f"Matched Exit Events in exit_data: {len(exit_data)}")
- print(f"Missing/NaN RR Exits: {exit_data['Marginal_RR'].isna().sum()}")
- # --- 4. Plotting Histograms per Round ---
- fig, axes = plt.subplots(3, 1, figsize=(10, 10), sharex=True, sharey=True)
- bins = np.linspace(-0.4, 0.4, 81)
- # Extract 3 shades of blue (Round 1: light blue, Round 2: medium blue, Round 3: dark blue)
- blues = sns.color_palette("Blues", n_colors=6)
- round_colors = [blues[2], blues[4], blues[5]]
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_exits = exit_data[exit_data['Round'] == roundN]['Marginal_RR'].dropna()
- if len(round_exits) > 0:
- mean_val = round_exits.mean()
- median_val = round_exits.median()
- sdDuration = round_exits.std()
- print(f'Round {roundN}:\nmean = {mean_val}\nmedian = {median_val}\nsd = {sdDuration}')
- sns.histplot(
- round_exits, bins=bins, kde=False, ax=ax,
- color=round_colors[i], alpha=0.85, stat="density",
- edgecolor="white", linewidth=0.8
- )
- sns.kdeplot(
- round_exits, fill=True,
- ax = ax, color=round_colors[i], alpha=0.45
- )
- # # Uncomment these if you want reference lines on the plot:
- # ax.axvline(0, color='gray', linestyle='--', linewidth=1.2, alpha=0.8, label='Equilibrium (0)')
- # ax.axvline(mean_val, color='#d62728', linestyle='-', linewidth=1.8, label=f'Mean ({mean_val:.3f})')
- # ax.axvline(median_val, color='#1f77b4', linestyle=':', linewidth=2.0, label=f'Median ({median_val:.3f})')
- ax.set_xlim(-0.4, 0.4)
- ax.set_title(f'Round {roundN} (N = {len(round_exits)} exits)', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Density', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3)
- sns.despine(ax=ax, top=True, right=True)
- axes[-1].set_xlabel('Marginal Reward Rate at Exit (reward/sec)', fontsize=11, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- # --- 1. Calculate Marginal RR ---
- resampled_df = resampledDF.copy()
- resampled_df['Marginal_RR'] = resampled_df['RewardRateCurrent'] - resampled_df['RewardRateTotal']
- exits_df = RealExitsDF.copy().rename(columns={'ParticipantID': 'ID'})
- # --- 2. Filter Out the First Exit of Each Round Per Participant ---
- # Sort to ensure chronological order per participant and round
- exits_df = exits_df.sort_values(['Round', 'ID', 'TimeStamp'])
- # Option A: If VisitNumber exists, keep visits greater than the minimum visit number
- if 'VisitNumber' in exits_df.columns:
- min_visits = exits_df.groupby(['Round', 'ID'])['VisitNumber'].transform('min')
- exits_df = exits_df[exits_df['VisitNumber'] > min_visits]
- # Option B: Fallback using row rank based on exit time
- else:
- exits_df['exit_rank'] = exits_df.groupby(['Round', 'ID'])['TimeStamp'].rank(method='first')
- exits_df = exits_df[exits_df['exit_rank'] > 1].drop(columns=['exit_rank'])
- # --- 3. Align Key Data Types Across DataFrames ---
- merge_by = ['Round', 'ID']
- if 'VisitNumber' in exits_df.columns and 'VisitNumber' in resampled_df.columns:
- merge_by.append('VisitNumber')
- for col in merge_by:
- exits_df[col] = exits_df[col].astype(int)
- resampled_df[col] = resampled_df[col].astype(int)
- # --- 4. Sort & Merge on Time Keys ---
- resampled_sorted = resampled_df.dropna(subset=['Marginal_RR']).sort_values('RegularTimestamp')
- exits_sorted = exits_df.sort_values('TimeStamp')
- exit_data = pd.merge_asof(
- exits_sorted,
- resampled_sorted[merge_by + ['RegularTimestamp', 'Marginal_RR']],
- left_on='TimeStamp', # Time column in RealExitsDF
- right_on='RegularTimestamp', # Time column in resampledDF
- by=merge_by,
- direction='nearest' # Matches to nearest floating timestamp
- )
- # Diagnostic verification
- print(f"Total Exit Events Filtered: {len(exits_df)}")
- print(f"Matched Exit Events in exit_data: {len(exit_data)}")
- print(f"Missing/NaN RR Exits: {exit_data['Marginal_RR'].isna().sum()}")
- # --- 5. Plotting Histograms per Round ---
- fig, axes = plt.subplots(3, 1, figsize=(10, 10), sharex=True, sharey=True)
- bins = np.linspace(-0.4, 0.4, 81)
- # Extract 3 shades of blue (Round 1: light blue, Round 2: medium blue, Round 3: dark blue)
- blues = sns.color_palette("Blues", n_colors=6)
- round_colors = [blues[2], blues[4], blues[5]]
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_exits = exit_data[exit_data['Round'] == roundN]['Marginal_RR'].dropna()
- if len(round_exits) > 0:
- mean_val = round_exits.mean()
- median_val = round_exits.median()
- sdDuration = round_exits.std()
- print(f'Round {roundN} (Excluding First Exits):\nmean = {mean_val}\nmedian = {median_val}\nsd = {sdDuration}\n')
- sns.histplot(
- round_exits, bins=bins, kde=False, ax=ax,
- color=round_colors[i], alpha=0.85, stat="density",
- edgecolor="white", linewidth=0.8
- )
- sns.kdeplot(
- round_exits, fill=True,
- ax=ax, color=round_colors[i], alpha=0.45
- )
- ax.set_xlim(-0.4, 0.4)
- ax.set_title(f'Round {roundN} (N = {len(round_exits)} exits - excluding first exit)', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Density', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3)
- sns.despine(ax=ax, top=True, right=True)
- axes[-1].set_xlabel('Marginal Reward Rate at Exit (reward/sec)', fontsize=11, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- # Filter resampledDF for nbhIn rows (where time in neighborhood starts)
- nbh_in_events = resampledDF[resampledDF['TimeInNbh'] < 2.0].sort_values('RegularTimestamp')
- exits_sorted = exit_data.sort_values('RegularTimestamp')
- # %%
- # --- 1. Standardize and Sort DataFrames ---
- nbh_df = nbh_in_events.copy()
- exits_df = RealExitsDF.copy().rename(columns={'ParticipantID': 'ID'})
- # Align data types
- for col in ['ID', 'Round']:
- nbh_df[col] = nbh_df[col].astype(int)
- exits_df[col] = exits_df[col].astype(int)
- # Sort chronologically within each participant & round
- nbh_sorted = nbh_df.sort_values(['ID', 'Round', 'RegularTimestamp']).reset_index(drop=True)
- exits_sorted = exits_df.sort_values(['ID', 'Round', 'TimeStamp']).reset_index(drop=True)
- # --- 2. Create Sequence Number per (ID, Round) Group ---
- nbh_sorted['visit_seq'] = nbh_sorted.groupby(['ID', 'Round']).cumcount()
- exits_sorted['visit_seq'] = exits_sorted.groupby(['ID', 'Round']).cumcount()
- # --- 3. Positional Inner Merge ---
- visit_durations_df = pd.merge(
- nbh_sorted,
- exits_sorted[['ID', 'Round', 'visit_seq', 'TimeStamp']].rename(columns={'TimeStamp': 'Exit_Timestamp'}),
- on=['ID', 'Round', 'visit_seq'],
- how='inner' # Drops trailing entries with no matching exit index
- )
- # --- 4. Calculate Duration & Ensure Exit > Entry ---
- visit_durations_df['VisitDuration'] = visit_durations_df['Exit_Timestamp'] - visit_durations_df['RegularTimestamp']
- # Filter out any non-positive duration edge cases
- visit_durations_df = visit_durations_df[visit_durations_df['VisitDuration'] > 0]
- # --- Verification Summary ---
- print(f"Total Entries: {len(nbh_in_events)}")
- print(f"Total Exits Available: {len(RealExitsDF)}")
- print(f"Successfully Matched Visits: {len(visit_durations_df)}")
- print(f"Dropped Uncompleted/Trailing Visits: {len(nbh_in_events) - len(visit_durations_df)}")
- # %%
- sns.set_theme(style="ticks")
- fig, axes = plt.subplots(3, 1, figsize=(10, 10), sharex=True, sharey=True)
- #bins = np.linspace(0, 200, 41)
- blues = sns.color_palette("Blues", n_colors=6)
- round_colors = [blues[2], blues[4], blues[5]] # [R1: Light, R2: Medium, R3: Dark]
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- round_durations = visit_durations_df[visit_durations_df['Round'] == roundN]['VisitDuration'].dropna()
- n_count = len(round_durations)
- meanDuration = np.mean(round_durations)
- medDuration = np.median(round_durations)
- sdDuration = np.std(round_durations)
- print(f'Round {roundN}:\nmean = {meanDuration}\nmedian = {medDuration}\nsd = {sdDuration}')
- sns.histplot(
- round_durations, kde=False, ax=ax, #bins=bins,
- color=round_colors[i], alpha=0.85, stat="density",
- edgecolor="white", linewidth=0.8
- )
- ax.set_title(f'Round {roundN} (N = {n_count} completed visits)', fontsize=12, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Density', fontsize=10, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3)
- sns.despine(ax=ax, top=True, right=True)
- axes[-1].set_xlabel(r'Neighborhood Visit Duration (seconds) [$t_{\mathrm{exit}} - t_{\mathrm{entry}}$]', fontsize=11, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- TIME_MAX = 400
- bins = np.linspace(0, TIME_MAX, 20)
- # --- 1. Blue Color Gradient across Rounds ---
- blues = sns.color_palette("Blues", n_colors=6)
- round_colors = [blues[2], blues[4], blues[5]]
- purps = sns.color_palette('Greys', n_colors=20)
- hist_colors = [purps[9], purps[14], purps[19]]
- # --- 2. Vectorized Data Processing ---
- plot_df = resampledDF.copy()
- # Find time of first valid payload per visit group
- valid_mask = plot_df['RewardRateCurrent'] > 0
- plot_df['first_payload_t'] = plot_df.where(valid_mask).groupby(['Round', 'ID', 'VisitNumber'])['TimeInNbh'].transform('min')
- # Mask timestamps prior to first payload (or visits with no payload)
- invalid_mask = (plot_df['TimeInNbh'] < plot_df['first_payload_t']) | plot_df['first_payload_t'].isna()
- plot_df.loc[invalid_mask, ['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- # Calculate Marginal Reward Rate (Patch RR - Expected RR) and filter binned time (<= 200s)
- plot_df['RR_Diff'] = plot_df['RewardRateCurrent'] - plot_df['RewardRateTotal']
- plot_df['TimeInNbh_bin'] = plot_df['TimeInNbh'].round(0)
- plot_df = plot_df[plot_df['TimeInNbh_bin'] <= TIME_MAX]
- # --- 3. Compact Visual Helpers ---
- def get_alphas(counts, min_val, max_val):
- """Rescale counts to alpha range [min_val, max_val]."""
- n_min, n_max = counts.min(), counts.max()
- if n_max > n_min:
- return min_val + (max_val - min_val) * (counts - n_min) / (n_max - n_min)
- return np.full_like(counts, max_val, dtype=float)
- def plot_fading_line(ax, x, y, counts, color, linewidth=2.5, zorder=2):
- alphas = get_alphas(counts, 0.20, 1.00)
- pts = np.array([x, y]).T.reshape(-1, 1, 2)
- segments = np.concatenate([pts[:-1], pts[1:]], axis=1)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- colors = [(*mcolors.to_rgb(color), a) for a in seg_alphas]
- ax.add_collection(LineCollection(segments, colors=colors, linewidths=linewidth, zorder=zorder))
- def plot_fading_ribbon(ax, x, y_lower, y_upper, counts, color, zorder=1):
- alphas = get_alphas(counts, 0.05, 0.25)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- base_rgb = mcolors.to_rgb(color)
- verts = [
- [[x[i], y_lower[i]], [x[i+1], y_lower[i+1]],
- [x[i+1], y_upper[i+1]], [x[i], y_upper[i]]]
- for i in range(len(x) - 1)
- ]
- ax.add_collection(PolyCollection(verts, facecolors=[(*base_rgb, a) for a in seg_alphas], edgecolors='none', zorder=zorder))
- # --- 4. Figure Plotting Loop ---
- fig, axes = plt.subplots(3, 1, figsize=(12, 20), sharex=True, sharey=True)
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- current_color = round_colors[i] # Progressive blue per round
- # Aggregate mean, sample size N, and compute 95% CI
- agg_df = (
- plot_df[plot_df['Round'] == roundN]
- .groupby('TimeInNbh_bin')['RR_Diff']
- .agg(mean='mean', sem='sem', count='count')
- .dropna(subset=['mean'])
- .reset_index()
- )
- # Calculate 95% Confidence Interval bound using scipy.stats.t
- t_crit = t_dist.ppf(0.975, df=np.maximum(agg_df['count'] - 1, 1))
- agg_df['ci95'] = (t_crit * agg_df['sem']).fillna(0)
- x = agg_df['TimeInNbh_bin'].values
- y_mean = agg_df['mean'].values
- ci = agg_df['ci95'].values
- counts = agg_df['count'].values
- # Plot elements
- ax.axhline(0, color='gray', linestyle='--', linewidth=1.2, alpha=0.7, zorder=0)
- plot_fading_ribbon(ax, x, y_mean - ci, y_mean + ci, counts, color=current_color, zorder=1)
- plot_fading_line(ax, x, y_mean, counts, color=current_color, linewidth=2.5, zorder=2)
- # Subplot Formatting
- ax.set_ylim(-0.05, 0.3)
- ax.set_xlim(0, TIME_MAX)
- ax.set_title(f'Round {roundN}', fontsize=24, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Marginal Reward Rate (reward/sec)', fontsize=20, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- axRight = ax.twinx()
- round_durations = visit_durations_df[visit_durations_df['Round'] == roundN]['VisitDuration'].dropna()
- n_count = len(round_durations)
- sns.histplot(
- round_durations, kde=False, ax=axRight, bins=bins,
- color=hist_colors[i], alpha=0.66, stat="density",
- edgecolor=None, linewidth=0.1
- )
- axRight.set_ylim(-.003, .018)
- axRight.set_ylabel('Density', fontsize=20, labelpad=8)
- axes[-1].set_xlabel('Time Since Neighborhood Entry (seconds)', fontsize=22, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- sns.set_theme(style="ticks")
- TIME_MAX = 400
- bins = np.linspace(0, TIME_MAX, 20)
- # --- 1. Blue Color Gradient across Rounds ---
- blues = sns.color_palette("Blues", n_colors=6)
- round_colors = [blues[2], blues[4], blues[5]]
- purps = sns.color_palette('Greys', n_colors=20)
- hist_colors = [purps[9], purps[14], purps[19]]
- # --- 2. Vectorized Data Processing ---
- plot_df = resampledDF.copy()
- # EXCLUDE FIRST PATCH VISIT (where current neighborhood entry == round's first entry)
- # plot_df = plot_df[plot_df['nbhIn'] != plot_df['firstIn']].copy()
- # Find time of first valid payload per visit group
- valid_mask = plot_df['RewardRateCurrent'] > 0
- plot_df['first_payload_t'] = plot_df.where(valid_mask).groupby(['Round', 'ID', 'VisitNumber'])['TimeInNbh'].transform('min')
- # Mask timestamps prior to first payload (or visits with no payload)
- invalid_mask = (plot_df['TimeInNbh'] < plot_df['first_payload_t']) | plot_df['first_payload_t'].isna()
- plot_df.loc[invalid_mask, ['RewardRateTotal', 'RewardRateCurrent']] = np.nan
- # Calculate Marginal Reward Rate (Patch RR - Expected RR) and filter binned time (<= 400s)
- plot_df['RR_Diff'] = plot_df['RewardRateCurrent'] - plot_df['RewardRateTotal']
- plot_df['TimeInNbh_bin'] = plot_df['TimeInNbh'].round(0)
- plot_df = plot_df[plot_df['TimeInNbh_bin'] <= TIME_MAX]
- plot_df = plot_df[plot_df['nbhIn'] != plot_df['firstIn']].copy()
- # --- 3. Compact Visual Helpers ---
- def get_alphas(counts, min_val, max_val):
- """Rescale counts to alpha range [min_val, max_val]."""
- n_min, n_max = counts.min(), counts.max()
- if n_max > n_min:
- return min_val + (max_val - min_val) * (counts - n_min) / (n_max - n_min)
- return np.full_like(counts, max_val, dtype=float)
- def plot_fading_line(ax, x, y, counts, color, linewidth=2.5, zorder=2):
- alphas = get_alphas(counts, 0.20, 1.00)
- pts = np.array([x, y]).T.reshape(-1, 1, 2)
- segments = np.concatenate([pts[:-1], pts[1:]], axis=1)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- colors = [(*mcolors.to_rgb(color), a) for a in seg_alphas]
- ax.add_collection(LineCollection(segments, colors=colors, linewidths=linewidth, zorder=zorder))
- def plot_fading_ribbon(ax, x, y_lower, y_upper, counts, color, zorder=1):
- alphas = get_alphas(counts, 0.05, 0.25)
- seg_alphas = 0.5 * (alphas[:-1] + alphas[1:])
- base_rgb = mcolors.to_rgb(color)
- verts = [
- [[x[i], y_lower[i]], [x[i+1], y_lower[i+1]],
- [x[i+1], y_upper[i+1]], [x[i], y_upper[i]]]
- for i in range(len(x) - 1)
- ]
- ax.add_collection(PolyCollection(verts, facecolors=[(*base_rgb, a) for a in seg_alphas], edgecolors='none', zorder=zorder))
- # --- 4. Figure Plotting Loop ---
- fig, axes = plt.subplots(3, 1, figsize=(12, 20), sharex=True, sharey=True)
- # Also exclude the first patch visit duration from the exits histogram dataset
- if 'firstIn' in visit_durations_df.columns and 'nbhIn' in visit_durations_df.columns:
- durations_filtered = visit_durations_df[visit_durations_df['nbhIn'] != visit_durations_df['firstIn']]
- elif 'VisitNumber' in visit_durations_df.columns:
- durations_filtered = visit_durations_df[visit_durations_df['VisitNumber'] > 1]
- else:
- durations_filtered = visit_durations_df.copy()
- for i, roundN in enumerate([1, 2, 3]):
- ax = axes[i]
- current_color = round_colors[i] # Progressive blue per round
- # Aggregate mean, sample size N, and compute 95% CI
- agg_df = (
- plot_df[plot_df['Round'] == roundN]
- .groupby('TimeInNbh_bin')['RR_Diff']
- .agg(mean='mean', sem='sem', count='count')
- .dropna(subset=['mean'])
- .reset_index()
- )
- # Calculate 95% Confidence Interval bound using scipy.stats.t
- t_crit = t_dist.ppf(0.975, df=np.maximum(agg_df['count'] - 1, 1))
- agg_df['ci95'] = (t_crit * agg_df['sem']).fillna(0)
- x = agg_df['TimeInNbh_bin'].values
- y_mean = agg_df['mean'].values
- ci = agg_df['ci95'].values
- counts = agg_df['count'].values
- # Plot elements
- ax.axhline(0, color='gray', linestyle='--', linewidth=1.2, alpha=0.7, zorder=0)
- plot_fading_ribbon(ax, x, y_mean - ci, y_mean + ci, counts, color=current_color, zorder=1)
- plot_fading_line(ax, x, y_mean, counts, color=current_color, linewidth=2.5, zorder=2)
- # Subplot Formatting
- ax.set_ylim(-0.05, 0.3)
- ax.set_xlim(0, TIME_MAX)
- ax.set_title(f'Round {roundN}', fontsize=24, fontweight='semibold', loc='left', pad=8)
- ax.set_ylabel('Marginal Reward Rate (reward/sec)', fontsize=20, labelpad=8)
- ax.grid(True, linestyle='--', alpha=0.3, zorder=0)
- sns.despine(ax=ax, top=True, right=True)
- axRight = ax.twinx()
- round_durations = durations_filtered[durations_filtered['Round'] == roundN]['VisitDuration'].dropna()
- sns.histplot(
- round_durations, kde=False, ax=axRight, bins=bins,
- color=hist_colors[i], alpha=0.66, stat="density",
- edgecolor=None, linewidth=0.1
- )
- axRight.set_ylim(-.003, .018)
- axRight.set_ylabel('Density', fontsize=20, labelpad=8)
- axes[-1].set_xlabel('Time Since Neighborhood Entry (seconds)', fontsize=22, labelpad=8)
- plt.tight_layout()
- plt.show()
- # %%
- AbandonHouseDF = VisitsDF[VisitsDF['AbandonedPatch'] == True]
- abandonIDs = AbandonHouseDF['ParticipantID'].unique()
- 1/30
- # %%
- AbandonParticipantsDF = VisitsDF[VisitsDF['ParticipantID'].isin(abandonIDs)]
- # %%
- R1PerfDF = PerformanceDF[PerformanceDF['Round'] == 1]
- R2PerfDF = PerformanceDF[PerformanceDF['Round'] == 2]
- R3PerfDF = PerformanceDF[PerformanceDF['Round'] == 3]
- # %%
- for df in [R1PerfDF, R2PerfDF, R3PerfDF]:
- print(df['PatchEfficiency'].mean())
- # %%
- def dobStringToAge(s):
- ## analysis date = 01/01/2025.
- ## returns approx. age as of this date.
- year = int(s[:4])
- month = int(s[4:6])
- day = int(s[6:])
- age = 2025 - (1 + year)
- extraDays = 360 - (30 * month) - day
- return age + extraDays / 360
- # %% [markdown]
- # ages = []
- # for dob in DemographicsDF['Date of Birth']:
- # dobTokens = dob.split('/')
- # dobString = ''
- # for t in dobTokens:
- # if len(t) == 1:
- # t = f'0{t}'
- # dobString = f'{t}{dobString}'
- # ages.append(dobStringToAge(dobString))
- #
- # #ages.append(dobStringToAge('20050424'))
- # #ages.append(dobStringToAge('20050803'))
- # print(sum(ages) / len(ages))
AnalysisControl.ipynb at commit 37e189e, no license · at the source
Overview
- Department of Psychology, Carnegie Mellon University, Pittsburgh, Pennsylvania, United States of America
- Laboratoire de Neurosciences Cognitives et Adaptatives (LNCA), Faculté de Psychologie, Université de Strasbourg, Strasbourg, France
- Entertainment Technology Center, Carnegie Mellon University, Pittsburgh, Pennsylvania, United States of America
- Neuroscience Institute, Carnegie Mellon University, Pittsburgh, Pennsylvania, United States of America
Abstract
The past decade has seen a resurgence of interest in the naturalistic behavior of patch foraging in both the neural and cognitive sciences. Patch foraging is observed throughout the animal kingdom and requires the integration of multiple cognitive processes, including value and perceptual learning, spatial estimation, the evaluation of instantaneous return rates, and anticipating future effort. Here, we present a novel virtual environment, BitterBuster, for studying patch foraging in humans. In the BitterBuster game, players navigate a virtual world with four neighborhoods (patches) populated with houses (harvest points) containing inhabitants known as “Oomplets”. The goal is to search for “bitter” Oomplets, characterized by the union of specific body colors and textures. Each neighborhood has a different overall probability of containing bitter Oomplets and, within each, specific house locations predict the relative rate of bitter Oomplets. BitterBuster can be used for solo or team based foraging, where the tasks of exploring and perceptual learning are implemented by two separate players working together. Our results show how both the solo and team version of BitterBuster can accurately measure multiple aspects of naturalistic patch foraging, with participants learning to identify optimal patches, classify harvestable resources, and make rational decisions about when to leave a patch before overharvesting. BitterBuster is well-suited for researchers interested in the contextual, neural, and cooperative dynamics of foraging.
Reproduced under the paper's license (CC BY), from the paper cited above.
Repositories
Its files are read in the Code ↔ Paper reader above, with 1 match between paragraphs and lines of code.
CoAxLab/BitterBusterTIER
37e189eddd719a3325713e82b5e5fb87134e3c28, 26 August 2026Availability: 1 check, the latest on 1 November 2026: the link answers
- 1 November 2026: the link answers
17 files
- Scripts/
AnalysisScripts/ , Jupyter, 421 lines.ipynb_checkpoints/ data_analysis_pipeline-c heckpoint.ipynb - Scripts/
AnalysisScripts/ , Jupyter, 2,002 lines, 1 matchAnalysisControl.ipynb - Scripts/
AnalysisScripts/ , Jupyter, 194 linesNavigation_Maps.ipynb - Scripts/
AnalysisScripts/ , Jupyter, 403 linesReward_vs_effort_clean.i pynb - Scripts/
AnalysisScripts/ , Jupyter, 1,118 linesdata_analysis_pipeline.i pynb - Scripts/
AnalysisScripts/ , Jupyter, 266 linessurvey_analysis.ipynb - Scripts/
ProcessingScripts/ , Python, 495 linesProcessFunctions.py - Scripts/
ProcessingScripts/ , Jupyter, 480 linesRawLogs_To_Analysis_tabl es.ipynb - Scripts/
ProcessingScripts/ , Python, 117 linesaddEffectiveReward.py - Scripts/
ProcessingScripts/ , Python, 23 linescreateDemographicsTable. py - Scripts/
ProcessingScripts/ , Python, 331 linescreatePerformanceTable.p y - Scripts/
ProcessingScripts/ , Python, 247 lineslogs_to_dfs.py - Scripts/
ProcessingScripts/ , Python, 186 linesprocess.py - Scripts/
ProcessingScripts/ , Jupyter, 299 linesprocessing.ipynb - Scripts/
ProcessingScripts/ , Python, 35 linesreplaceIDs.py - Scripts/
ProcessingScripts/ , Jupyter, 60 linessurvey_process_spreadshe et.ipynb - README.md, Text, 13 lines
CoAxLab/BitterBusterForagingPaperRelease
e95090abb8147060239c9ac5ad53b2237f790206, 20 June 2024Availability: 1 check, the latest on 1 November 2026: the link answers
- 1 November 2026: the link answers
2 files
- LICENSE.md, License, 173 lines
- README.md, Text, 9 lines
The paper's code and data availability statement is in the Data section.
Tracing map
Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.
What the map holds:
- 2 repositories of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
- 16 scripts, each with its path and the digest of its content;
- 1 match between paragraphs of the paper and lines of the code (method lexical-v1);
- neither the text of the paper nor the code itself.
Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.
Data
No dataset and no data link were found in the paper.
Data Availability
The data set and accompanying code is available at GitHub via https://
Reproduced under the paper's license (CC BY), from the paper cited above.
Versions
The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.
Version 1, 1 November 2026: the first record
Recorded: type, language, journal, volume, issue, pages, dates, 8 authors, 5 MeSH terms, 1 funder, 19 references.
Cite
This paper
Kasarda, J. P., Bahuguna, J., Zhang, A., Tong, H., Tan, Y., Wang, R., Tarr, M. J., & Verstynen, T. (2026). BitterBuster: A virtual world for understanding naturalistic patch foraging. PloS one, 21(10), e0361502. https://
BibTeX
@article{kasarda2026bitt
author = {Kasarda, John P. and Bahuguna, Jyotika and Zhang, Angela and Tong, Hua and Tan, Yuan and Wang, Ruizi and Tarr, Michael J. and Verstynen, Timothy},
title = {{BitterBuster: A virtual world for understanding naturalistic patch foraging}},
journal = {PloS one},
year = {2026},
month = oct,
volume = {21},
number = {10},
pages = {e0361502},
publisher = {PLOS},
issn = {1932-6203},
doi = {10.1371/
url = {https://
pmid = {42961977},
pmcid = {PMC13722190}
}
RIS
TY - JOUR
AU - Kasarda, John P.
AU - Bahuguna, Jyotika
AU - Zhang, Angela
AU - Tong, Hua
AU - Tan, Yuan
AU - Wang, Ruizi
AU - Tarr, Michael J.
AU - Verstynen, Timothy
TI - BitterBuster: A virtual world for understanding naturalistic patch foraging
T2 - PloS one
J2 - PLoS One
PY - 2026
DA - 2026/
VL - 21
IS - 10
SP - e0361502
SN - 1932-6203
PB - PLOS
DO - 10.1371/
UR - https://
LA - en
ER -
CSL-JSON
{
"id": "10.1371/
"type": "article-journal",
"title": "BitterBuster: A virtual world for understanding naturalistic patch foraging",
"container-title": "PloS one",
"author": [
{
"family": "Kasarda",
"given": "John P."
},
{
"family": "Bahuguna",
"given": "Jyotika"
},
{
"family": "Zhang",
"given": "Angela"
},
{
"family": "Tong",
"given": "Hua"
},
{
"family": "Tan",
"given": "Yuan"
},
{
"family": "Wang",
"given": "Ruizi"
},
{
"family": "Tarr",
"given": "Michael J."
},
{
"family": "Verstynen",
"given": "Timothy"
}
],
"container-title-short":
"volume": "21",
"issue": "10",
"page": "e0361502",
"DOI": "10.1371/
"PMID": "42961977",
"PMCID": "PMC13722190",
"ISSN": "1932-6203",
"publisher": "PLOS",
"URL": "https://
"language": "en",
"issued": {
"date-parts": [
[
2026,
10,
30
]
]
}
}
The tracing map gets a citation of its own once an author has validated it and it has a DOI.
Similar papers
The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.
- [1] doi:10.1038/s41467-026-68711-x [code]
- Neural and computational mechanisms underlying one-shot perceptual learning in humans.Journal: Nature communicationsIn common: statannotations, NetworkX, Pillow, 7 other tools, cognitive
- [2] doi:10.1016/j.isci.2026.115488 [code]
- An integrated &
lt;i& gt;i& lt;/ i& gt; & lt;i& gt;n vitro& lt;/ i& gt; platform and biophysical modeling approach for studying synaptic transmission in isolated neuronal pairs. Journal: iScienceIn common: statannotations, NetworkX, Pillow, 7 other tools - [3] doi:10.1038/s41586-025-09520-y [code]
- ABCA7 variants impact phosphatidylcholine and mitochondria in neurons.Journal: NatureIn common: statannotations, NetworkX, Pillow, 7 other tools
- [4] doi:10.1093/schizbullopen/sgae010 [code]
- Cortical Network Disruption Is Minimal in Early Stages of Psychosis.Journal: Schizophrenia bulletin openIn common: statannotations, NetworkX, Pillow, 7 other tools
- [5] doi:10.1038/s41598-023-40683-8 [code]
- Cross-species comparative analysis of single presynapses.Journal: Scientific reportsIn common: statannotations, NetworkX, Pillow, 7 other tools
- [6] doi:10.1371/journal.pcbi.1012629 [code]
- Image-based force inference by biomechanical simulation.Journal: PLoS computational biologyIn common: statannotations, NetworkX, Pillow, 7 other tools
- [7] doi:10.1371/journal.pcbi.1012989 [code]
- Foraging animals use dynamic Bayesian updating to model meta-uncertainty in environment representations.Journal: PLoS computational biologyIn common: NetworkX, seaborn, pandas, 3 other tools, 3 references
- [8] doi:10.1038/s42003-024-06858-3 [code]
- Dynamical structure-function correlations provide robust and generalizable signatures of consciousness in humans.Journal: Communications biologyIn common: statannotations, NetworkX, statsmodels, 6 other tools, cognitive
- [9] doi:10.7554/elife.101688 [code]
- Transcranial direct current stimulation modulates primate brain dynamics across states of consciousness.Journal: eLifeIn common: statannotations, NetworkX, statsmodels, 6 other tools, cognitive
- [10] doi:10.3389/fnagi.2025.1550673 [code]
- Disrupted calcium dynamics and electrophysiological activity in the stratum pyramidale and hippocampal alveus during fear conditioning in the 5xFAD model of Alzheimer's disease.Journal: Frontiers in aging neuroscienceIn common: statannotations, NetworkX, statsmodels, 6 other tools, cognitive
Contribute
The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.
Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.
Claim this paper
Correct its record
Say what each link of this record is, remove the ones that are not the paper's, add the ones that are missing. The correction becomes a new version of the record, in its Versions section.
Validate its tracing map
You validate the map as this page shows it: 2 repositories of the authors' code, each at its verified commit and with its license, 16 scripts, and 1 match between paragraphs and code (see the Code and Map sections). It then receives a DOI on Zenodo, with you (your ORCID iD) and OSCR as its creators; the code itself is not deposited.
The map's fingerprint: sha256:41894a7c9df74220…
Add the badge to the code's README
If you are an author or a maintainer of this paper's code, you can add this badge. It says what OSCR knows of the paper (its code, its validated tracing map) and links the code to this page. Paste one of these snippets under the README's title. Keep the two marks around it, so that the badge can be replaced in place and is never added twice. OSCR changes nothing in your repository.
Markdown
<!-- oscr-badge:start -->
[](https://openscicode.org/badge/doi/10.1371/journal.pone.0361502)
<!-- oscr-badge:end -->reStructuredText
.. oscr-badge:start
.. image:: https://openscicode.org/badge/doi/10.1371/journal.pone.0361502.svg
:target: https://openscicode.org/badge/doi/10.1371/journal.pone.0361502
:alt: OSCR
.. oscr-badge:endHTML
<!-- oscr-badge:start -->
<a href="https://openscicode.org/badge/doi/10.1371/journal.pone.0361502"><img src="https://openscicode.org/badge/doi/10.1371/journal.pone.0361502.svg" alt="OSCR"></a>
<!-- oscr-badge:end -->On GitHub, open the README, choose Edit (the pencil), paste the snippet under the title, then "Commit changes".
- CoAxLab/BitterBusterTIER: README.md on GitHub
- CoAxLab/BitterBusterForagingPaperRelease: README.md on GitHub
If you choose to, you can have the badge kept up to date for you, with oscr badge --write from the command line or with the GitHub Action in your repository's workflows. Every other form and option (a card, AsciiDoc, LaTeX, shields.io) is on the badge page.
Request its removal
To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).
Discussion, reproductions, activity
Discussion: code errors and mismatches between the code and the paper about this paper and its code, from readers and its authors, beside the paper and its tracing map.
Report a mismatch between the code and the paper · Report a code error (signed in)
Reproductions: reports from readers who ran the authors' code and did not get the paper's result, what they ran, with which environment, commit and data.
Report a reproduction failure (signed in)
Versions of its code: the releases of its code tied to a version of this paper (the preprint, the accepted manuscript, the version of record…), each with the tracing map versioned with it.
Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.
