finish main scripts and some unit tests

Author: mdhumphries

Functions/Summaries_of_Beta_distribution.py

@@ -4,7 +4,7 @@
 
 def summaries_of_Beta_Distribution(Alpha, Beta, stat_type, *args):
     if np.isnan(Alpha) or np.isnan(Beta):
-        statistic = None
+        statistic = np.nan
     else:
         if stat_type == 'MAP':
             x = np.arange(0, 1, 0.001)

Functions/__pycache__/Summaries_of_Beta_distribution.cpython-39.pyc

@@ -0,0 +1,36 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Wed Aug 24 12:14:59 2022
+
+@author: Mark Humphries
+"""
+import numpy as np
+from Functions.Summaries_of_Beta_distribution import summaries_of_Beta_Distribution
+
+# trial_data is a dict of {Alpha,Beta,MAPprobability,Precision}
+# dataframe_row_of_prior_trials is a Pandas dataframe row of the previous trial
+# (or the empty dataframe on the first trial)
+
+def interpolate_null_trials(dict_of_trial_data,dataframe_row_of_prior_trial,alpha_zero,beta_zero):
+    
+    # if not NAN, then no need to interpolate: just assign the supdated alpha etc to the interpolated variables
+    if not np.isnan(dict_of_trial_data['Alpha']):
+        interpolated_values = {'Alpha_interpolated': dict_of_trial_data['Alpha'], 'Beta_interpolated': dict_of_trial_data['Beta'], \
+                      'MAPprobability_interpolated': dict_of_trial_data['MAPprobability'], 'Precision_interpolated': dict_of_trial_data['Precision']}
+       
+      
+    elif dataframe_row_of_prior_trial.empty:
+    # it is a null trial: check if dataframe row is empty - then is first trial, so assign priors
+         MAPprobability = summaries_of_Beta_Distribution(alpha_zero,beta_zero,'MAP')
+         precision = summaries_of_Beta_Distribution(alpha_zero,beta_zero,'precision') 
+         interpolated_values = {'Alpha_interpolated': alpha_zero, 'Beta_interpolated': beta_zero, \
+                       'MAPprobability_interpolated': MAPprobability, 'Precision_interpolated': precision}
+    else:
+    # it is a null trial, and the dataframe is not empty, so assign interpolated values of previous trial
+          interpolated_values = {'Alpha_interpolated': dataframe_row_of_prior_trial.at['Alpha_interpolated'], 'Beta_interpolated': dataframe_row_of_prior_trial.at['Beta_interpolated'], \
+                        'MAPprobability_interpolated': dataframe_row_of_prior_trial.at['MAPprobability_interpolated'], 'Precision_interpolated': dataframe_row_of_prior_trial.at['Precision_interpolated']}  
+    
+    
+    dict_of_trial_data.update(interpolated_values)  # add interpolated values to others from this trial
+    return dict_of_trial_data # a dict of all results for this trial

Functions/__pycache__/interpolate_null_trials.cpython-39.pyc

@@ -1,3 +1,4 @@
+import numpy as np
 
 def update_strategy_posterior_probability(trial_type, decay_rate, success_total, failure_total, alpha0, beta0):
     if trial_type == "success":
@@ -11,8 +12,8 @@ def update_strategy_posterior_probability(trial_type, decay_rate, success_total,
         alpha = (alpha0 + success_total)
         beta = (beta0 + failure_total)
     else:
-        alpha = None
-        beta = None
+        alpha = np.nan
+        beta = np.nan
 
     return success_total, failure_total, alpha, beta
 

Functions/__pycache__/update_strategy_posterior_probability.cpython-39.pyc

@@ -8,12 +8,13 @@
 from Functions.update_strategy_posterior_probability import update_strategy_posterior_probability
 from Functions.Summaries_of_Beta_distribution import summaries_of_Beta_Distribution
 from Functions.plotSessionStructure import plotSessionStructure
+from Functions.interpolate_null_trials import interpolate_null_trials
 
 # initiate TestData variable so that rat 2 testdata can be loaded
 TestData = pd.read_csv('data.csv')
 
 #  choose strategies to evaluate: subset shown in Figure 1
-strategies = ['go_left','go_right']
+strategies = ['go_left','go_right','go_cued','win_stay_spatial','lose_shift_cued','lose_shift_spatial']
 
 # set prior
 prior_type = 'Uniform' #set prior type
@@ -27,7 +28,7 @@
 event_totals = {}   # empty dict to store totals of events for each strategy
 # initialise dataframes
 for index_strategy in range(len(strategies)):
-    Output_collection[strategies[index_strategy]] =  pd.DataFrame(columns = ['Alpha', 'Beta', 'MAPprobability', 'Precision'])  # empty Dataframe to input data into
+    Output_collection[strategies[index_strategy]] =  pd.DataFrame(columns = ['Alpha', 'Beta','MAPprobability', 'Precision','Alpha_interpolated', 'Beta_interpolated','MAPprobability_interpolated', 'Precision_interpolated'])  # empty Dataframe to input data into
     event_totals[strategies[index_strategy]] = {}; # create empty dict for this strategy
     event_totals[strategies[index_strategy]]['success_total'] = 0;
     event_totals[strategies[index_strategy]]['failure_total'] = 0;
@@ -49,25 +50,53 @@
 
        MAPprobability = summaries_of_Beta_Distribution(Alpha, Beta, 'MAP')
        precision = summaries_of_Beta_Distribution(Alpha, Beta, 'precision')
-    
+      
+              
+       #%% interpolate null trials
+       this_trials_data= {'Alpha':Alpha, 'Beta':Beta, 'MAPprobability':MAPprobability, 'Precision':precision}     # create dict of current data to pass
+       if trial > 0:
+           previous_trials_data = Output_collection[strategies[index_strategy]].iloc[trial-1] 
+       else:
+           previous_trials_data = Output_collection[strategies[index_strategy]]  # pass empty dataframe
+           
+       new_row_of_data = interpolate_null_trials(this_trials_data,previous_trials_data,alpha0,beta0)
+       
        # store results  - dynamically-defined dataframe...
-       new_row = {'Alpha':Alpha, 'Beta':Beta, 'MAPprobability':MAPprobability, 'Precision':precision}     # create new row for dataframe as a dict
-       new_df= pd.DataFrame([new_row])   # have to convert to dataframe to use concat!!
+       new_df= pd.DataFrame([new_row_of_data])   # have to convert to dataframe to use concat!!
        Output_collection[strategies[index_strategy]] = pd.concat([Output_collection[strategies[index_strategy]], new_df], ignore_index=True)       # add new row to dataframe
-       
+      
+      
 #%% plot results
 no_Trials = np.size(TestData.TrialIndex)
 
-# plotting time series of MAPprobability
+# plotting time series of MAPprobability for Rule Strategies
 plt.figure(figsize=(10, 5))
 plt.plot(Output_collection['go_left'].MAPprobability, linewidth=0.75)  # plots the time series
-plt.axis([0, no_Trials, 0, 1.25])  # establishes axis limits
+plt.plot(Output_collection['go_right'].MAPprobability, linewidth=0.75, color=(0.4, 0.8, 0.5))  # plots the time series
+plt.plot(Output_collection['go_cued'].MAPprobability, linewidth=0.75, color=(0.8,0.6,0.5))  # plots the time series
+plt.axis([0, no_Trials, 0, 1.25])  # establishes axis limits 
 plt.xlabel('Trials'), plt.ylabel('P(Strategy)')  # labelling the axis
-plt.axhline(y=0.5, color='firebrick', linestyle='--', linewidth=0.75,
-            label="Chance")  # shows the line at which Chance is exceeded
+plt.axhline(y=0.5, color='darkgrey', linewidth=0.75, label="Chance")  # shows the line at which Chance is exceeded
 
 plotSessionStructure(TestData)
-
 plt.show()       
+
+# plotting Precision for the same three strategies (precision identical for go_left and go_right)
+plt.figure(figsize=(10, 5))
+plt.plot(Output_collection['go_left'].Precision, linewidth=0.75)  # plots the time series
+plt.plot(Output_collection['go_cued'].Precision, linewidth=0.75, color=(0.8,0.6,0.5))  # plots the time series
+plt.xlabel('Trials'), plt.ylabel('Precision')  # labelling the axis
+plt.show()      
 
-       
+
+# plotting MAP probability for some exploratory strategies - use interpolated values
+plt.figure(figsize=(10, 5))
+plt.plot(Output_collection['lose_shift_cued'].MAPprobability_interpolated, linewidth=0.75, color=(1, 0.1, 0.6))  # plots the time series
+plt.plot(Output_collection['lose_shift_spatial'].MAPprobability_interpolated, linewidth=0.75, color=(0.8, 0.6, 0.5))  # plots the time series
+plt.plot(Output_collection['win_stay_spatial'].MAPprobability_interpolated, linewidth=0.75, color=(0.4,0.8,0.5))  # plots the time series
+plt.axis([0, no_Trials, 0, 1.25])  # establishes axis limits 
+plt.xlabel('Trials'), plt.ylabel('P(Strategy)')  # labelling the axis
+plt.axhline(y=0.5, color='darkgrey', linewidth=0.75, label="Chance")  # shows the line at which Chance is exceeded
+
+plotSessionStructure(TestData)
+plt.show()  

Functions/interpolate_null_trials.py

@@ -84,7 +84,7 @@ def alternate(rows):
     # checks if the subject made a different choice on this trial from the previous obe
     nTrials = len(rows)
      # "at" selects the value at the row/column location in the dataframe
-    if nTrials > 1 & rows.at[nTrials-1,'Choice'] != rows.at[nTrials-2,'Choice']:      # check the current trial's choice
+    if nTrials > 1 and rows.at[nTrials-1,'Choice'] != rows.at[nTrials-2,'Choice']:      # check the current trial's choice
         trial_type = "success"
     else:
         trial_type = "failure"
@@ -99,11 +99,11 @@ def lose_shift_cued(rows):
     nTrials = len(rows)
      # "at" selects the value at the row/column location in the dataframe
      # check that the previous trial was not rewarded ('lose')
-    if nTrials > 1 & rows.at[nTrials-2,'Reward'] == "no":   
+    if nTrials > 1 and rows.at[nTrials-2,'Reward'] == "no":   
         # now check if the subject shifted their cued-based choice
-        if rows.at[nTrials-2,'Choice'] == rows.at[nTrials-2,'CuePosition'] & rows.at[nTrials-1,'Choice'] != rows.at[nTrials-1,'CuePosition']:
+        if rows.at[nTrials-2,'Choice'] == rows.at[nTrials-2,'CuePosition'] and rows.at[nTrials-1,'Choice'] != rows.at[nTrials-1,'CuePosition']:
             trial_type = "success"  # shifted from cued to uncued choice
-        elif rows.at[nTrials-2,'Choice'] != rows.at[nTrials-2,'CuePosition'] & rows.at[nTrials-1,'Choice'] == rows.at[nTrials-1,'CuePosition']:
+        elif rows.at[nTrials-2,'Choice'] != rows.at[nTrials-2,'CuePosition'] and rows.at[nTrials-1,'Choice'] == rows.at[nTrials-1,'CuePosition']:
             trial_type = "success" # shifted from uncued to cued chpice
         else:
             trial_type = "failure"
@@ -119,7 +119,7 @@ def lose_shift_spatial(rows):
     nTrials = len(rows)
      # "at" selects the value at the row/column location in the dataframe
      # check that the previous trial was not rewarded ('lose')
-    if nTrials > 1 & rows.at[nTrials-2,'Reward'] == "no":   
+    if nTrials > 1 and rows.at[nTrials-2,'Reward'] == "no":   
         # now check if the subject shifted their spatial choice
         if rows.at[nTrials-1,'Choice'] != rows.at[nTrials-2,'Choice']: 
             trial_type = "success"
@@ -134,7 +134,7 @@ def sticky(rows):
     # checks if the subject made the same choice on this trial as the previous one
     nTrials = len(rows)
      # "at" selects the value at the row/column location in the dataframe
-    if nTrials > 1 & rows.at[nTrials-1,'Choice'] == rows.at[nTrials-2,'Choice']:      # check the current trial's choice
+    if nTrials > 1 and rows.at[nTrials-1,'Choice'] == rows.at[nTrials-2,'Choice']:      # check the current trial's choice
         trial_type = "success"
     else:
         trial_type = "failure"
@@ -149,11 +149,11 @@ def win_stay_cued(rows):
     nTrials = len(rows)
      # "at" selects the value at the row/column location in the dataframe
      # check that the previous trial was rewarded ('win')
-    if nTrials > 1 & rows.at[nTrials-2,'Reward'] == "yes":   
+    if nTrials > 1 and rows.at[nTrials-2,'Reward'] == "yes":   
         # now check if the subject stayed with the cued-based choice
-        if rows.at[nTrials-2,'Choice'] == rows.at[nTrials-2,'CuePosition'] & rows.at[nTrials-1,'Choice'] == rows.at[nTrials-1,'CuePosition']:
+        if rows.at[nTrials-2,'Choice'] == rows.at[nTrials-2,'CuePosition'] and rows.at[nTrials-1,'Choice'] == rows.at[nTrials-1,'CuePosition']:
             trial_type = "success"  # made the same cued choice
-        elif rows.at[nTrials-2,'Choice'] != rows.at[nTrials-2,'CuePosition'] & rows.at[nTrials-1,'Choice'] != rows.at[nTrials-1,'CuePosition']:
+        elif rows.at[nTrials-2,'Choice'] != rows.at[nTrials-2,'CuePosition'] and rows.at[nTrials-1,'Choice'] != rows.at[nTrials-1,'CuePosition']:
             trial_type = "success" # made the same uncued choice
         else:
             trial_type = "failure"
@@ -169,7 +169,7 @@ def win_stay_spatial(rows):
     nTrials = len(rows)
      # "at" selects the value at the row/column location in the dataframe
      # check that the previous trial was rewarded ('win')
-    if nTrials > 1 & rows.at[nTrials-2,'Reward'] == "yes":   
+    if nTrials > 1 and rows.at[nTrials-2,'Reward'] == "yes":   
         # now check if the subject stayed with the same spatial choice
         if rows.at[nTrials-1,'Choice'] == rows.at[nTrials-2,'Choice']: 
             trial_type = "success"

Functions/update_strategy_posterior_probability.py

@@ -0,0 +1,69 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Wed Aug 24 19:39:45 2022
+
+Unit tests for the strategy models
+
+@author: Mark Humphries
+"""
+
+import pandas as pd
+import strategymodels
+
+# make test dataframe
+test_data_for_rules = {'Choice' : ['left','left','left','left','right','right','right','right'],
+        'CuePosition': ['left','right','left','right','left','right','left','right'],
+        'Reward': ['yes','yes','no','no','yes','yes','no','no']}
+TestData = pd.DataFrame(test_data_for_rules)
+
+# target results for rules
+
+go_left_target = ['success','success','success','success','failure','failure','failure','failure']
+go_right_target = ['failure','failure','failure','failure','success','success','success','success']
+go_cued_target = ['success','failure','success','failure','failure','success','failure','success']
+go_uncued_target = ['failure','success','failure','success','success','failure','success','failure']
+
+#%% run rule strategy models
+trial_type_left = []; trial_type_right = []; trial_type_cued = []; trial_type_uncued = [];
+for trial in range(len(TestData)):
+    rows_of_data = TestData.iloc[0:trial+1]     # select all rows of data up to the curren trial; is trial+1 as dataframe includes column row as row 0????
+    trial_type_left.append(strategymodels.go_left(rows_of_data))         # test whether go-left was used       
+    trial_type_right.append(strategymodels.go_right(rows_of_data))
+    trial_type_cued.append(strategymodels.go_cued(rows_of_data))  
+    trial_type_uncued.append(strategymodels.go_uncued(rows_of_data))  
+    
+#%% did they pass the tests?    
+print('Go left passed?')    
+print(trial_type_left == go_left_target)
+
+print('Go right passed?')    
+print(trial_type_right == go_right_target)
+
+print('Go cued passed?')    
+print(trial_type_cued == go_cued_target)
+    
+print('Go uncued passed?')    
+print(trial_type_uncued == go_uncued_target)
+
+#%%
+#  test Explore strategy models
+# 
+
+# extend test data to include all win-stay and lose-shift test cases
+test_data_for_rules['Choice'].append('right','left','right','right','left')
+test_data_for_rules['CuePosition'].append('right','right','left','left','left')
+test_data_for_rules['Reward'].append('yes','yes','yes','no','no')
+
+TestData = pd.DataFrame(test_data_for_rules)  # overwrite dataframe
+
+# define target results for exploration rules (null, success, failure)
+
+
+# run explore strategy models
+trial_type_alternate = []; trial_type_sticky = []; trial_type_win_stay_spatial= []; 
+trial_type_win_stay_cue = []; trial_type_lose_shift_spatial= []; trial_type_lose_shift_cue = []; 
+
+for trial in range(len(TestData)):
+    rows_of_data = TestData.iloc[0:trial+1]     # select all rows of data up to the curren trial; is trial+1 as dataframe includes column row as row 0????  
+