The following code has been used to compare couples of vectors in a given structure composed of 4×12 cells which result from simulated results of 4×12 model genotypes (download archive here). Most vectors have a similar dimension (175), but some differ, requiring a multiple random comparison by permutation. The final distribution of the p-values (reported in the four subplots) provides an indication of likelihood of statistical significant difference.
The first part of the code loads the data structure, then generates two random vectors (randperm) that are used to create the permutation of the target data to compare (which are selected assigning a value between 1 and 12 to the parameters el1 and el2). The operation is repeated n time (default: 500) to allow multiple (n) permutations.
The second part of the code is used to provide instructions for a series of charts. This is meant to show how to change some basic settings such as color of the histograms, axis range, histogram edge etc. By assigning the output of the function histogram to a variable (default: fig), it is possible to display the content of this variable and access the function names assigned by Matlab which are used to specify further plotting settings. E.g.:
fig = 1
Histogram with properties:Data: [500×1 double]
Values: [3 3 10 22 27 34 27 30 26 28 38 26 38 28 29 30 33 41 27]
NumBins: 19
BinEdges: [1×20 double]
BinWidth: 0.0500
BinLimits: [0.0500 1]
Normalization: ‘count’
FaceColor: [0.8000 0.2000 0.6000]
EdgeColor: [0.8000 0.2000 0.2000]
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 |
% TTest comparison % NB Check for multiple comparisons: e.g. Bonferroni correction % https://en.wikipedia.org/wiki/Bonferroni_correction clear all close all load('DrugGoalRatioBPf2 Ratio in_904-999 End Drug Assumption First Period.mat') el1=7; %first element in comparison el2=8; %second element in comparison for i=1:500 for k=1:4 n=length(ratio{k,el1}); %lengtn of first vector to compare m=length(ratio{k,el2}); %lengtn of second vector to compare poll=min([n m]); %length of vetor to generate permutation of indexes select1=randperm(n, poll); %index permutation vec 1 select2=randperm(m, poll); %index permutation vec 2 vec1 =ratio{k,el1}(select1); vec2 =ratio{k,el2}(select2); [X P CI STATS]=ttest(vec1, vec2); store_b(i,k)=P; % if k==1 % disp(mean(vec1)); % disp(mean(vec2)); % end end end % percentage of p<0.05 passb=store_b<0.05; for k=1:4 percb(k)=100*sum(passb(:,k))/length(passb(:,k)); end % generates a character array with the found P values p_titles=char( num2str(percb(1)), num2str(percb(2)), num2str(percb(3)), num2str(percb(4))); %plotting with extra settings h=figure; screen_size(3:4)=[1800 900]; set(h,'Position',screen_size) subplot(2,2,1) fig=histogram(store_b(:,1), 'BinEdges', 0:0.05:1); set(gca,'XTick',0:.05:1, 'FontWeight','bold', 'FontSize',10) set(gca,'YTick',0:10:1500, 'YTickLabel',0:1:200, 'FontWeight','bold', 'FontSize',10) topb=max(fig.Values); axis([0 1 0 topb+(topb/5)]) text(0.05,topb+(topb/10),[' percP = % ',p_titles(1,:)],... 'VerticalAlignment','middle',... 'HorizontalAlignment','left',... 'FontSize',16) subplot(2,2,2) fig=histogram(store_b(:,2), 'BinWidth', 0.05); set(gca,'XTick',0:.05:1, 'FontWeight','bold', 'FontSize',10) set(gca,'YTick',0:10:1500, 'YTickLabel',0:1:200, 'FontWeight','bold', 'FontSize',10) topb=max(fig.Values); axis([0 1 0 topb+(topb/5)]) text(0.05,topb+(topb/10),[' percP = % ',p_titles(2,:)],... 'VerticalAlignment','middle',... 'HorizontalAlignment','left',... 'FontSize',16) subplot(2,2,3) fig=histogram(store_b(:,3), 'BinWidth', 0.05, 'FaceColor', 'g'); set(gca,'XTick',0:.05:1, 'FontWeight','bold', 'FontSize',10); set(gca,'YTick',0:10:1500, 'YTickLabel',0:1:200, 'FontWeight','bold', 'FontSize',10); topb=max(fig.Values); axis([0 1 0 topb+(topb/5)]) text(0.05,topb+(topb/10),[' percP = % ',p_titles(3,:)],... 'VerticalAlignment','middle',... 'HorizontalAlignment','left',... 'FontWeight', 'bold', ... 'Color', 'r', ... 'FontSize',16) subplot(2,2,4) fig=histogram(store_b(:,3), 'BinWidth', 0.05, 'FaceColor', [0.8 0.2 0.6], 'EdgeColor', [0.8 0.2 0.2]); set(gca,'XTick',0:.05:1, 'FontWeight','bold', 'FontSize',10); set(gca,'YTick',0:10:1500, 'YTickLabel',0:1:200, 'FontWeight','bold', 'FontSize',10); topb=max(fig.Values); axis([0 1 0 topb+(topb/5)]) text(0.05,topb+(topb/10),[' percP = % ',p_titles(4,:)],... 'VerticalAlignment','middle',... 'HorizontalAlignment','left',... 'FontWeight', 'bold', ... 'Color', [0.2, 0.2, 0.8], ... 'FontSize',16) |