Chi square test for fraction of runs/sessions per animal. USe python script...

Chi square test for fraction of runs/sessions per animal. USe python script since it does not qutomatically merge groups

matlab/+Script/+data/number_of_betas_per_animal.m

@@ -8,11 +8,65 @@ beta_tables_taste = util.load.all_taste_betas().combined;
 %% return locations of which set for all rows in s
 X = beta_tables_taste(:,{'only_animal_idx', 'tastant'});
 s = X.Variables;
-[~,ia,ic]=unique(s,'rows'); 
+[~,ia,ic]=unique(s,'rows');
 n = histc(ic,unique(ic));
 T = table(s(ia,1),s(ia,2),n);
 % unstack the tabe
 T_unstacked = unstack(T,'n','Var2');
-T_unstacked = renamevars(T_unstacked, 'Var1', 'animal_idx');
+T_unstacked.Properties.RowNames = cellstr(T_unstacked.Var1);
+T_unstacked.Var1 = [];
+
+% convert nans to zeros
+var_vals = T_unstacked.Variables;
+var_vals(isnan(var_vals)) = 0;
+T_unstacked.Variables = var_vals;
+% sort animals
+T_unstacked = sortrows(T_unstacked,1,'descend');
+
+%% Stacked bar plot
+fh = figure; ax = gca;
+fractions = (T_unstacked.Variables./sum(T_unstacked.Variables,2))*100;
+ba = bar(fractions,...
+    'stacked','FaceColor','flat');
+colors = util.plot.getColors().taste;
+
+for i=1:length(ba)
+    ba(i).CData = colors(i,:);
+end
+xticklabels(T_unstacked.Properties.RowNames)
+
+util.plot.cosmeticsSave...
+    (fh, ax, 3.2, 1.4, cur_save_dir,...
+     'barplot_fractions.svg', 'off','off');
 %% Save the table
-save(fullfile(cur_save_dir, 'betas_taste_animal.mat'), 'T_unstacked');
\ No newline at end of file
+save_mat = false;
+if save_mat
+    save(fullfile(cur_save_dir, 'betas_taste_animal.mat'), 'T_unstacked');
+end
+
+%% IGNORE(USE PYTHON chi_square_taste_fractions.py): Chi square proportions test
+% first calculated expected numbers for each animal based on the sum of the
+% session numbers
+total_proportions = sum(T_unstacked.Variables,1)./...
+    sum(T_unstacked.Variables,'all');
+observed_per_animal = T_unstacked.Variables;
+expected_per_animal = sum(T_unstacked.Variables,2) * total_proportions;
+num_animals = length(observed_per_animal);
+
+for i = 1:num_animals
+    [h(i),p(i),stats] = chi2gof([1,2,3],'Frequency',observed_per_animal(i, :),...
+        'Expected',expected_per_animal(i, :),'Alpha',0.05/num_animals);
+    disp(stats)
+end
+%% Write the csv for python
+% NOTE: matlab automatically pools the data from multiple taste groups when
+% the n is small for the last three animals, use python which does not pool
+% automatically (chi_square_taste_fractions.py)
+writetable(array2table(observed_per_animal,...
+    'VariableNames',T_unstacked.Properties.VariableNames,...
+    'RowNames',T_unstacked.Properties.RowNames),...
+    fullfile(cur_save_dir, 'observed.csv'),'WriteRowNames',true)
+writetable(array2table(expected_per_animal,...
+    'VariableNames',T_unstacked.Properties.VariableNames,...
+    'RowNames',T_unstacked.Properties.RowNames),...
+    fullfile(cur_save_dir, 'expected.csv'),'WriteRowNames',true)

matlab/+util/+plot/getColors.m

@@ -16,5 +16,7 @@ colors.blue = [46,49,144]./255;
 colors.modules = [0.941176470588235,0.501960784313726,0.501960784313726;...
     0.250980392156863,0.878431372549020,0.815686274509804;...
     0.254901960784314,0.411764705882353,0.882352941176471];
+
+colors.taste = [201, 145, 192; 123, 173, 213; 212, 219, 149]./255;
 end
 

python/brainNetwork/scripts/chi_square_taste_fractions.py

+from numpy import array_split
+from scipy.stats import chisquare
+import pandas as pd
+from pathlib import Path
+# read observed and expected frequency per taste
+main_dir = '/Users/karimia/code/Brain_Network/matlab/figs/betas_taste_animal'
+dirs = [Path(main_dir).joinpath(x+'.csv') for x in ['observed', 'expected']]
+Ts = [pd.read_csv(dir,index_col=0) for dir in dirs]
+arrays = [t.to_numpy() for t in Ts]
+# Chi square test
+chiq, p_vals = chisquare(f_obs=arrays[0], f_exp=arrays[1], axis=1)
+significance = p_vals < (0.05/8)
+
+# Create table and display
+p_value_T = pd.DataFrame([chiq, p_vals, significance],columns=Ts[0].index,
+                         index=['chi_test_statistic', 'p_vals_chi_square', 'significance_bonf']).T
+p_value_T.to_csv(Path(main_dir,'chi_square_animal_taste.csv'))
\ No newline at end of file