123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101 |
- % Appies HSMM-MVPA as in Anderson 2016 to find the boundaries of cognitive
- % processing stages
- bumpM = 4; % number of bumps, 4 bumps == 5 stages
- nBuMax = max(bumpM);
- nModel = length(bumpM);
- nPca = 32; % 32 PCA -> more than 95% variance, 25 -> more than 85%; 20 -> 80% variance on PCA
- nSj = 18;
- rep = 150;
- shape = 2; % gamma distributoin, shape parameter
- roiDataR = 'nScoreAveR'; % 'nScoreAve' or 'nScorePca'
- dataPath = './data';
- dataRead = '/vars4hsmmBL_R.mat';
- roiTypeR = ['roi' roiDataR(end-3:end)];
- roiTypeR(end) = [];
- saveName = ['/hmmBu' num2str(bumpM) 'R' roiTypeR 'Map3StagBL.mat'];
- addpath('./hsmm/');
- info = {['test ', num2str(bumpM),' bumps n repetitions one models: R; map 4 stage, PCA cov by subject, Baseline']};
- load([dataPath, dataRead], roiDataR, 'condR','xR','yR','subjtR')
- % right hand
- data_ = eval(roiDataR);
- data = single(data_(:,1:nPca));
- x = xR;
- y = yR;
- condR_ = condR;
- condR_(logical( condR==2 | condR==3 )) = 1; % fan1-target
- condR_(logical( condR==4 | condR==5 )) = 2; % fan2-target
- condR_(logical( condR==6 | condR==7 )) = 3; % fan1-foil
- condR_(logical( condR==8 | condR==9 )) = 4; % fan2-foil
- cond = condR_;
- subj = subjtR;
- % left hand
- clear data_ roiDataR xR yR condR condR_ subjtR
- % define map settings, 3th stage (retrieval) different gamma
- map = [1, 1, 1, 1, 1;
- 1, 1, 2, 1, 1;
- 1, 1, 3, 1, 1;
- 1, 1, 4, 1, 1];
- gamIni = repmat([2 200], nBuMax+1,1);
- magIni = -1 + 2*rand(nPca, nBuMax, rep);
- % bump kernel correlation, removed form inside HSMM
- data = calcBumps(data{1,1});
- lkhr = zeros(repnModel); % [repetitions, hand response, models(number bumps),]
- parpool(30);
- parfor r = 1:rep
- % do HSMM
- for m = 1:nModel % iterate models
- [lkhr(r,m),~,~,~]=hsmmEEG50CondMapAllDNew(data,cond,magIni(:,1:bumpM(m),r),gamIni(1:bumpM(m)+1,:),1,x,y,map,shape);
- end
- end
- end
- disp('End of random initializations')
- % recompute best models
- lkh = cell(1,nModel); % hand, model(nBump)
- bumpMag = cell(1,nModel);
- flatGam = cell(1,nModel);
- proBump = cell(1,nModel);
- for m = 1:nModel
- [~,iBest] = max(lkhr(:,m));
- [lkh{1,m},bumpMag{1,m},flatGam{1,m},proBump{1,m}] = hsmmEEG50CondMapAllDNew(data,cond,magIni(:,1:bumpM(m),iBest),gamIni(1:bumpM(m)+1,:),1,x,y,map,shape);
- end
- % LOOCV, leave-one-out cross validation
- lkhCV = zeros(nSj,nModel);
- parfor sj = 1:nSj
- for m = 1:nModel % iterate models
- sTrain = logical(subj ~= sj);
- sTest = logical(subj == sj);
- [~, magTrain, gammTrain, ~] = hsmmEEG50CondMapAllDNew(data,cond(sTrain),bumpMag{1,m},flatGam{1,m}, 1 ,x{h}(sTrain),y(sTrain),map,shape);
- [ lkhCV(sj, m), ~, ~, ~] = hsmmEEG50CondMapAllDNew(data,cond(sTest), magTrain, gammTrain, 0 ,x(sTest),y(sTest),map,shape);
- end
- end
- %save(['./dataSave' saveName], 'lkhr','lkh','bumpMag','flatGam','proBump','info','nPca', 'rep', 'x','y','cond','lkhCV','subj')
- %format longG
- disp(char(datetime(now,'ConvertFrom','datenum')))
- delete(gcp('nocreate'))
- clear all
- close all
|