private/calculatePST.m (90f4202) - SVMCrossVal.git

calculatePST.m

function pst = calculatePST(des,globalStart,baselineStart,baselineEnd,globalEnd,eventList,data,sessionList)
    bstart          = baselineStart;
    bend            = baselineEnd;
    edur            = 12;
    pre             =  globalStart;
    post            =  globalEnd;
    res             = 1;

normz           = 'file';
    pm              = 0;

lsess           = getNumberOfScans(des);
    nSessions       = getNumberOfSessions(des);
    tr              = 2;

[evntrow evntcol]=size(eventList);

hsec=str2num(des.xsDes.High_pass_Filter(8:end-3)); % Highpass filter [sec] from SPM.mat

if strcmp(des.xBF.UNITS,'secs')
        unitsecs=1;
    end;

nScansPerSession=getNumberOfScans(des);
    %stime=[0:tr:max(nScansPerSession)*tr+post-tr]; % Stimulus time for raw data plot
    stime=0:tr:max(nScansPerSession)*tr+round(post/tr)*tr-tr; % Stimulus time for raw data plot

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % RUN
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Digital Highpass
    Rp=0.5;
    Rs=20;
    NO=1;
    Wp=1/((1/2/tr)/(1/hsec));
    [B, A] = ellip(NO,Rp,Rs,Wp,'high');

sdata(1:max(nScansPerSession)+round(post/tr),1:nSessions)=nan; % Open Data Matrix
    for z=1:nSessions % Fill Data Matrix sessionwise
        sdata(1:nScansPerSession(z),z)=data(sum(nScansPerSession(1:z))-nScansPerSession(z)+1:sum(nScansPerSession(1:z)))';
    end;
%         usdata=sdata; % Keep unfiltered data

sdatamean=nanmean(nanmean(sdata(:,:)));
    for z=1:nSessions
%             X(:,z)=[1:1:max(nScansPerSession)]'; % #Volume
        sdata(1:nScansPerSession(z),z)=filtfilt(B,A,sdata(1:nScansPerSession(z),z)); %Filter Data (Highpass)
    end;
    sdata=sdata+sdatamean;

%%%%Parametric Modulation Modus%%%%
    if pm %Find Parameters for Event of Interest
        [imods modss mods erow evntrow eventList] = getParametricMappingEvents(eventList,evntrow,des,pmf);
    end;
    %%%%PM%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

for zr=1:evntrow
        n{zr}=0;
        nn{zr}=0; 
        nnn{zr}=0;
        sstart{zr}=1;
    end;

sesst0=0; 
    for sessionID=sessionList
        if sessionID>1
            sesst0(sessionID)=sum(lsess(1:sessionID-1))*tr;  
        end;
        for zr=1:evntrow  %LABEL NUMBER, EVENT GROUP
            sstart{zr}=n{zr}+1;
            for ze=1:evntcol % EVENT INDEX in EventList
                if ze==1 || (ze>1 && eventList(zr,ze)~=eventList(zr,ze-1))
                    for zz=1:length(des.Sess(sessionID).U(eventList(zr,ze)).ons) % EVENT REPETITION NUMBER
                        if ~unitsecs
                            des.Sess(sessionID).U(eventList(zr,ze)).ons(zz)=(des.Sess(sessionID).U(eventList(zr,ze)).ons(zz)-1)*tr;
                            des.Sess(sessionID).U(eventList(zr,ze)).dur(zz)=(des.Sess(sessionID).U(eventList(zr,ze)).dur(zz)-1)*tr;
                        end;

nnn{zr}=nnn{zr}+1; % INFO for rawdataplot start
                        if des.Sess(sessionID).U(eventList(zr,ze)).dur(zz)<edur
                            mev{zr}(nnn{zr},1:2)=[des.Sess(sessionID).U(eventList(zr,ze)).ons(zz)+sesst0(sessionID) edur]; % modeled event [onset length]
                        else
                            mev{zr}(nnn{zr},1:2)=[des.Sess(sessionID).U(eventList(zr,ze)).ons(zz)+sesst0(sessionID) des.Sess(sessionID).U(eventList(zr,ze)).dur(zz)];
                        end; % INFO for rawdataplot end

n{zr}=n{zr}+1;
                        pst{zr}(n{zr},:)=interp1(stime,sdata(:,sessionID),[des.Sess(sessionID).U(eventList(zr,ze)).ons(zz)+pre:res:des.Sess(sessionID).U(eventList(zr,ze)).ons(zz)+post],'linear');
                        if strcmp(normz,'epoc')
                            bline=nanmean(pst{zr}(n{zr},round(-pre/res+(bstart)/res+1):round(-pre/res+(bend)/res+1)));
                            if isnan(bline)
                                pst{zr}(n{zr},1:-pre/res+post/res+1)=nan;
                            else
%                                     nn{zr}=nn{zr}+1;
                                pst{zr}(n{zr},:)=(pst{zr}(n{zr},:)-bline)/bline*100; % 'epoch-based' normalization
                            end;
                        end;
                    end;
                end;
            end;
            if ~strcmp(normz,'epoc')
                bline(zr)=nanmean(nanmean(pst{zr}(sstart{zr}:n{zr},-pre/res+(bstart)/res+1:-pre/res+(bend)/res+1)));
                bstd(zr)=nanmean(nanstd(pst{zr}(sstart{zr}:n{zr},-pre/res+(bstart)/res+1:-pre/res+(bend)/res+1)));
                nn{zr}=n{zr};
            end;
        end;
        if strcmp(normz,'filz')
            for zr=1:evntrow
                pst{zr}(sstart{zr}:n{zr},:)=(pst{zr}(sstart{zr}:n{zr},:)-mean(bline))/mean(bstd); % session-based z-score normalization
            end;
        elseif strcmp(normz,'file')
            for zr=1:evntrow
                pst{zr}(sstart{zr}:n{zr},:)=(pst{zr}(sstart{zr}:n{zr},:)-mean(bline))/mean(bline)*100; % session-based normalization
            end;
        end;
    end;
    
    %%%%%%%%%%% new 090109 Axel: "Normalization" for SVM
    
%     for zr=1:evntrow
% %        tmp_max(zr)=max(nanmean(pst{zr}));
%           tmp_max(zr)=max(max(pst{zr}));
%           tmp_min(zr)=min(min(pst{zr}));
%     end;
% %    pstmax=max(tmp_max);
% %    pstmin=min(tmp_min);
%     pstmax=max(abs([tmp_max tmp_min]));
%     for zr=1:evntrow
%         pst{zr}=pst{zr}./pstmax-.5;
% %        pst{zr}=(pst{zr}-pstmin)./(pstmax-pstmin);
%     end;
    
    
    
    
end

SVMCrossVal.git

plot optimiert. normalisierung pst, zwischenversion UI Defaults verbessert