|
| 1 | + |
| 2 | +function acrOut = aggregateAcr(acrIn, varargin) |
| 3 | +% function acrOut = aggregateAcr(acrIn, varargin) |
| 4 | +% |
| 5 | +% Translate voxel acronyms into a useful list to split analysis by. |
| 6 | +% |
| 7 | +% e.g.: |
| 8 | +% >> acrIn = {'VISp5', 'VISp6a', 'MB', 'VISpm1', 'MB', 'LP', 'VISp5', 'LGd', 'VISp'}; |
| 9 | +% >> acrOut = aggregateAcr(acrIn); acrOut |
| 10 | +% acrOut = |
| 11 | +% 1×9 cell array |
| 12 | +% {'VISp'} {'VISp'} {[NaN]} {'VISpm'} {[NaN]} {'LP'} {'VISp'} {'LGd'} {'VISp'} |
| 13 | + |
| 14 | +% Optional input argument aggType can be 1 for a short list of high-level |
| 15 | +% regions, or 2 [default] for a long list of sensible regions (manually |
| 16 | +% selected by N.A.S.). E.g. |
| 17 | +% >> acrOut = aggregateAcr(acrIn, 1); acrOut |
| 18 | +% acrOut = |
| 19 | +% 1×9 cell array |
| 20 | +% {'Isocortex'} {'Isocortex'} {'MB'} {'Isocortex'} {'MB'} {'TH'} {'Isocortex'} {'TH'} {'Isocortex'} |
| 21 | + |
| 22 | +% This function will return NaN for the acronyms that don't have a parent |
| 23 | +% in the list or that aren't valid acronyms. To find these in the output, |
| 24 | +% use: |
| 25 | +% >> cellfun(@(x)any(isnan(x)), acrOut) |
| 26 | + |
| 27 | +% If you are calling this function many times, you can speed up execution |
| 28 | +% by providing the tree-pair representation of the structure tree, which is |
| 29 | +% a bit slow to calculate. Provide this as the third input argument. You |
| 30 | +% can make it like this: |
| 31 | +% >> [~, tp] = makeSTtree(st); |
| 32 | + |
| 33 | +st = loadStructureTree(); |
| 34 | + |
| 35 | +if nargin>2 |
| 36 | + tp = varargin{2}; |
| 37 | +else |
| 38 | + % see doc of this function to know what tp is |
| 39 | + [~, tp] = makeSTtree(st); |
| 40 | +end |
| 41 | + |
| 42 | +if nargin>1 |
| 43 | + aggType = varargin{1}; |
| 44 | + if isempty(aggType); aggType = 2; end |
| 45 | +else |
| 46 | + aggType = 2; |
| 47 | +end |
| 48 | + |
| 49 | +aggList = getAggList(aggType); |
| 50 | + |
| 51 | +aggTP = tp(ismember(tp(:,1), aggList),:); |
| 52 | +% aggTP now contains only children of desired list |
| 53 | +% since these lists are chosen such that every entry of structure tree only |
| 54 | +% has one (at most) parent on this list, the second column now contains |
| 55 | +% just one entry for each of our desired areas. |
| 56 | + |
| 57 | +% we also want exact matches so let's just add those in here. |
| 58 | +aggTP = [aggTP; horzcat(aggList', aggList', zeros(size(aggList))')]; |
| 59 | + |
| 60 | + |
| 61 | +if ~iscell(acrIn) % can be either a single acronym string, or a cell array |
| 62 | + acrIn = {acrIn}; |
| 63 | +end |
| 64 | + |
| 65 | +allAcr = st.acronym; |
| 66 | +allId = st.id; |
| 67 | + |
| 68 | +notAnAcr = ~ismember(acrIn, allAcr); |
| 69 | +if any(notAnAcr) |
| 70 | + firstBadAcr = acrIn{find(notAnAcr,1)}; |
| 71 | + warning(... |
| 72 | + sprintf('One or more of your acronyms didn''t parse. One that didn''t is ''%s''.', ... |
| 73 | + firstBadAcr)); |
| 74 | +end |
| 75 | + |
| 76 | +acrInClean = acrIn(~notAnAcr); % only process the good ones, the rest we set to NaN later |
| 77 | + |
| 78 | +idIn = cellfun(@(x)allId(strcmp(allAcr,x)), acrInClean); |
| 79 | + |
| 80 | +tpLoc = arrayfun(@(x)find(aggTP(:,2)==x), idIn, 'uni', false); |
| 81 | + |
| 82 | +noMatch = cellfun(@isempty, tpLoc); |
| 83 | +multipleMatch = cellfun(@(x)numel(x)>1, tpLoc); |
| 84 | + |
| 85 | +if any(multipleMatch) |
| 86 | + % -- todo -- report one that is bad |
| 87 | + warning('Some acronyms appear to be children of two or more parents in the list. This probably means your aggList is badly constructed. Will return the first match.'); |
| 88 | + tpLoc = cellfun(@(x)x(1), tpLoc, 'uni', false); |
| 89 | +end |
| 90 | + |
| 91 | +% not throwing a warning here because this is probably common/normal, e.g. |
| 92 | +% voxels that are "MB" or "root" or other things. |
| 93 | +% These are returned as NaN. |
| 94 | +tpLocClean = tpLoc(~noMatch); |
| 95 | + |
| 96 | +tpLocClean = cell2mat(tpLocClean); |
| 97 | + |
| 98 | +aggParent = aggTP(tpLocClean,1); |
| 99 | + |
| 100 | +aggOut = arrayfun(@(x)allAcr{allId==x}, aggParent, 'uni', false); |
| 101 | + |
| 102 | +% now put in NaN's for the ones that didn't work |
| 103 | +acrOut1 = cell(size(tpLoc)); |
| 104 | +acrOut1(noMatch) = {NaN}; |
| 105 | +acrOut1(~noMatch) = aggOut; |
| 106 | + |
| 107 | +acrOut = cell(size(acrIn)); |
| 108 | +acrOut(notAnAcr) = {NaN}; |
| 109 | +acrOut(~notAnAcr) = acrOut1; |
| 110 | + |
| 111 | +function aggList = getAggList(aggType) |
| 112 | + |
| 113 | +switch aggType |
| 114 | + case 1 % 10 top-level divisions: {'Isocortex', 'HPF', 'OLF', 'CTXsp', 'CNU', 'TH', 'HY', 'MB', 'HB', 'CB'} |
| 115 | + aggList = [315 1089 698 703 623 549 1097 313 1065 512]; |
| 116 | + |
| 117 | + case 2 % 316 manually selected regions |
| 118 | + aggList = [184 985 993 353 329 337 345 369 361 182305689 ... |
| 119 | + 378 1057 677 1011 1002 1027 1018 402 394 ... |
| 120 | + 409 385 425 533 312782574 312782628 39 48 972 44 723 ... |
| 121 | + 731 738 746 104 111 119 894 879 886 312782546 417 ... |
| 122 | + 541 922 895 507 151 159 597 605 814 961 619 639 647 ... |
| 123 | + 788 566 382 423 463 726 982 19 918 926 843 1037 1084 ... |
| 124 | + 502 484682470 589508447 484682508 583 952 966 131 295 ... |
| 125 | + 319 780 672 56 998 754 250 258 266 310 333 23 292 536 ... |
| 126 | + 1105 403 1022 1031 342 298 564 596 581 351 629 685 718 ... |
| 127 | + 725 733 741 563807435 406 609 1044 475 170 218 1020 1029 ... |
| 128 | + 325 560581551 255 127 64 1120 1113 155 59 362 366 1077 ... |
| 129 | + 149 15 181 560581559 189 599 907 575 930 560581563 262 ... |
| 130 | + 27 563807439 178 321 483 186 390 38 30 118 ... |
| 131 | + 223 72 263 272 830 452 523 1109 126 133 347 286 338 ... |
| 132 | + 576073699 689 88 210 491 525 557 515 980 1004 63 693 ... |
| 133 | + 946 194 226 364 576073704 173 470 614 797 302 4 580 271 ... |
| 134 | + 874 381 749 607344830 246 128 294 795 215 531 ... |
| 135 | + 628 634 706 1061 549009203 616 214 35 549009211 975 115 ... |
| 136 | + 606826663 757 231 66 75 58 374 1052 12 100 197 591 872 ... |
| 137 | + 612 7 867 398 280 880 599626927 898 931 1093 318 534 574 ... |
| 138 | + 621 549009215 549009219 549009223 549009227 679 147 162 ... |
| 139 | + 604 146 238 350 358 207 96 101 711 1039 903 642 651 429 ... |
| 140 | + 437 445 589508451 653 661 135 839 1048 372 83 136 106 ... |
| 141 | + 203 235 307 395 852 859 938 177 169 995 1069 209 202 225 ... |
| 142 | + 217 765 773 781 206 230 222 912 976 984 1091 936 944 951 ... |
| 143 | + 957 968 1007 1056 1064 1025 1033 1041 1049 989 91 846 ... |
| 144 | + 589508455]; |
| 145 | + |
| 146 | +end |
0 commit comments