diff --git a/usedMcode/callAllPos.m b/usedMcode/callAllPos.m index 7eb61a1..6cf43e6 100644 --- a/usedMcode/callAllPos.m +++ b/usedMcode/callAllPos.m @@ -1,25 +1,28 @@ addpath('/home/hohlochj/masterarbeit/usedMcode') pathToFile='/nfs/wsi/ti/messor/hohlochj/origData/'; -%pathToFile='/home/jph/Uni/masterarbeit/Block1_ReachingMovements/'; +%pathToFile='/home/jph/Uni/masterarbeit/origData/'; maxFile=5; threshold=10000; %EMG is classified as movement EMGChannels={'AbdPolLo','Biceps','Triceps','FrontDelt','MidDelt','BackDelt'}; +noSynergies=3; -allSubjects=true; %run all subjects and days or only one random day for one random subject - +allSubjects=false; %run all subjects and days or only one random day for one random subject +name='defaultPos'; %suffix for the output, has to be valid name for file (no space, /, ...) windowEMG=0.2; windowEEG=1; -shift=0.2; -pburgOrder=50; -minEEGFreq=0; -maxEEGFreq=200; +shiftEMG=0.05; +shiftEEG=0.2; +eegOffset=0; %predict actions x*shiftEEG after EEG measurement +pburgOrder=250; +minEEGFreq=2; +maxEEGFreq=49; pause=false; -ridgeParams=[0.1]; +noLFsamples=5; +ridgeParams=100; k=10; maxExpC=0; maxPerClass=250; -EEG=true; -poolObj=parpool(22); +poolObj=parpool(3); [subjects,numbers]=namesAndNumbers(pathToFile); numbersMat=cell2mat(numbers); @@ -27,51 +30,46 @@ j=1; for i=1:size(subjects,2) subject=subjects{i}; - for number=numbers{i} + for number=numbers{i} subjectsForNumbers{j}=subject; j=j+1; end end -j=j-1; +j=j-1; %number of trial-days for all subjects if allSubjects + maxRidgeParamIndexEEGpos=zeros([j,3,k]); + maxRidgeParamIndexEMGpos=zeros([j,3,k]); + maxRidgeParamIndexLFpos=zeros([j,3,k]); + correlationEMGpos=zeros([j,3]); %x,y,angle + correlationEEGpos=zeros([j,3]); + correlationLFpos=zeros([j,3]); - %meanAccurancysEMG=zeros([j,1]); - %meanAccurancysEEG=zeros([j,1]); - %maxCEMG=zeros([j,k,1]); - %maxCEEG=zeros([j,k,1]); - cmScaledEMG=zeros([j,5,5]); - cmScaledEEG=zeros([j,5,5]); - maxRidgeParamIndex=zeros([j,3,k]); - correlationEMG=zeros([j,3]); %x,y,angle - correlationEEG=zeros([j,3]); - parfor j=1:size(numbersMat,2) number=numbersMat(j); subject=subjectsForNumbers{j}; - savePath=readAllPos(pathToFile,subject,number,windowEMG,windowEEG,shift,maxFile,threshold,pburgOrder,minEEGFreq,maxEEGFreq,pause,EMGChannels); - % [meanAccurancysEMG(j),maxCEMG(j,:),cmScaledEMG(j,:,:)]=svmEciton(savePath,~EEG,k,maxExpC,maxPerClass); - % [meanAccurancysEEG(j),maxCEEG(j,:),cmScaledEEG(j,:,:)]=svmEciton(savePath,EEG,k,maxExpC,maxPerClass); - [correlationEMG(j,:),maxRidgeParamIndex(j,:,:)]=ridgeCV(savePath,false,k,ridgeParams); - [correlationEEG(j,:),maxRidgeParamIndex(j,:,:)]=ridgeCV(savePath,EEG,k,ridgeParams); - fprintf('%s%i finished %s\n',subject,number,datestr(datetime('now'))) + savePath=readAllPos(pathToFile, subject,number,windowEMG,windowEEG,shiftEMG,shiftEEG,maxFile,... + threshold,pburgOrder,minEEGFreq,maxEEGFreq,pause,noLFsamples,EMGChannels); + [correlationEMGpos(j,:),maxRidgeParamIndex(j,:,:)]=ridgeCV(savePath,'EMG','kin',k,ridgeParams,eegOffset); + [correlationEEGpos(j,:),maxRidgeParamIndex(j,:,:)]=ridgeCV(savePath,'EEG','kin',k,ridgeParams,eegOffset); + [correlationLFpos(j,:),maxRidgeParamIndex(j,:,:)]=ridgeCV(savePath,'LF','kin',k,ridgeParams,eegOffset); + %fprintf('%s%i finished %s\n',subject,number,datestr(datetime('now'))) end - - save(strcat(pathToFile,sprintf('../matlabData/%s_callAllPos.mat',datestr(datetime('now'))))); + + save(strcat(pathToFile,sprintf('../matlabData/%s_callAll-%s.mat',datestr(datetime('now')),name))); else j=fix(rand()*size(numbersMat,2)+1); number=numbersMat(j); subject=subjectsForNumbers{j}; - savePath=readEEG(pathToFile,subject,number,windowEMG,windowEEG,shift,maxFile,threshold,pburgOrder,minEEGFreq,maxEEGFreq,pause); - %[meanAccurancysEMG,maxCEMG,cmScaledEMG(:,:)]=svmEciton(savePath,~EEG,k,maxExpC,maxPerClass); - %[meanAccurancysEEG,maxCEEG,cmScaledEEG(:,:)]=svmEciton(savePath,EEG,k,maxExpC,maxPerClass); - [correlationEMG,maxRidgeParamIndex(:,:)]=ridgeCV(savePath,false,k,ridgeParams); - [correlationEEG,maxRidgeParamIndex(:,:)]=ridgeCV(savePath,EEG,k,ridgeParams); - fprintf('%s%i finished %s\n',subject,number,datestr(datetime('now'))) - - save(strcat(pathToFile,sprintf('../matlabData/%s_call%s%iPos.mat',datestr(datetime('now')),subject,number))); + savePath=readAllPos(pathToFile, subject,number,windowEMG,windowEEG,shiftEMG,shiftEEG,maxFile,... + threshold,pburgOrder,minEEGFreq,maxEEGFreq,pause,noLFsamples,EMGChannels); + [correlationEMGpos,maxRidgeParamIndex(:,:)]=ridgeCV(savePath,'EMG','kin',k,ridgeParams,eegOffset); + [correlationEEGpos,maxRidgeParamIndex(:,:)]=ridgeCV(savePath,'EEG','kin',k,ridgeParams,eegOffset); + [correlationLFpos,maxRidgeParamIndex(:,:)]=ridgeCV(savePath,'LF','kin',k,ridgeParams,eegOffset); + + save(strcat(pathToFile,sprintf('../matlabData/%s_call%s%i-%s.mat',datestr(datetime('now')),subject,number,name))); end delete(poolObj) diff --git a/usedMcode/generateTrainingDataPos.m b/usedMcode/generateTrainingDataPos.m index 597ce31..6bb7d97 100644 --- a/usedMcode/generateTrainingDataPos.m +++ b/usedMcode/generateTrainingDataPos.m @@ -1,34 +1,27 @@ -function [trainingDataEEG,trainingDataEEGlf,trainingDataEMG,kinPerSec] = generateTrainingDataPos(signal,kin,windowEMG,windowEEG,shift,params,pburgOrder,minEEGFreq,maxEEGFreq,noLFsamples,EMGChannels) +function [trainingDataEEG,trainingDataEEGlf,trainingDataEMG,kinPerSec] = generateTrainingDataPos(signal,kin,windowEMG,windowEEG,shiftEMG,shiftEEG,params,pburgOrder,minEEGFreq,maxEEGFreq,noLFsamples,EMGChannels) bci_sf=params.SamplingRate.NumericValue; - signalWindowEMG=bci_sf*windowEMG; signalWindowEEG=bci_sf*windowEEG; - tempEEG=zeros([32,fix(floor(size(signal,1)/signalWindowEEG-1)*windowEEG/shift+1),maxEEGFreq-minEEGFreq+1]); - tempEEGlf=zeros([32,fix(floor(size(signal,1)/signalWindowEEG-1)*windowEEG/shift+1),noLFsamples]); - tempEMG=zeros(fix([6,fix(size(signal,1)/signalWindowEMG-1)*windowEMG/shift+1])); + tempEEG=zeros([32,fix(floor(size(signal,1)/signalWindowEEG-1)*windowEEG/shiftEEG+1),201]); + tempEEGlf=zeros([32,fix(floor(size(signal,1)/signalWindowEEG-1)*windowEEG/shiftEEG+1),noLFsamples]); %Filter around 50Hz and below 2 Hz [A,B]= butter(2,[48 52]/(bci_sf/2),'stop'); - [C,D]= butter(2,1/(bci_sf/2),'high'); + [C,D]= butter(2,minEEGFreq/(bci_sf/2),'high'); [E,F]= butter(2,[148 152]/(bci_sf/2),'stop'); + [G,H]= butter(2,maxEEGFreq/(bci_sf/2),'low'); % filter for low Frequencies [V,W]= butter(2,0.01/(bci_sf/2),'high'); [X,Y]= butter(2,1/(bci_sf/2),'low'); parfor i=1:32 %filter single channel, w/o EMG, HEOG, Synchro and 0s - tempEEG(i,:,:)=shiftingPburg(filtfilt(double(E),double(F),filtfilt(double(C),double(D),filtfilt(double(A),double(B),double(signal(:,i))))),bci_sf,windowEEG,shift,pburgOrder,minEEGFreq,maxEEGFreq,EMGChannels); - tempEEGlf(i,:,:)=shiftingDownsample(filtfilt(double(X),double(Y),filtfilt(double(V),double(W),double(signal(:,i)))),bci_sf,windowEEG,shift,noLFsamples); + tempEEG(i,:,:)=shiftingPburg(filtfilt(double(G),double(H),filtfilt(double(E),double(F),filtfilt(double(C),double(D),filtfilt(double(A),double(B),double(signal(:,i)))))),bci_sf,windowEEG,shiftEEG,pburgOrder,0,200); + tempEEGlf(i,:,:)=shiftingDownsample(filtfilt(double(X),double(Y),filtfilt(double(V),double(W),double(signal(:,i)))),bci_sf,windowEEG,shiftEEG,noLFsamples); end trainingDataEEG=permute(tempEEG,[2 1 3]); trainingDataEEGlf=permute(tempEEGlf,[2 1 3]); - usedChannels={'AbdPolLo','Biceps','Triceps','FrontDelt','MidDelt','BackDelt'}; - emgchannels=find(ismember(params.ChannelNames.Value,usedChannels)); - - parfor i=1:6; - tempEMG(i,:)=waveformLength(filtfilt(double(E),double(F),filtfilt(double(C),double(D),filtfilt(double(A),double(B),double(signal(:,emgchannels(i)))))),bci_sf,windowEMG,shift); - end - trainingDataEMG=permute(tempEMG,[2 1 3]); + trainingDataEMG=waveformLengthAll(signal(:,ismember(params.ChannelNames.Value,EMGChannels)),bci_sf,windowEMG,shiftEMG); %shift kin according to synch channel in first minute [~,startingPoint]=max(abs(diff(signal(1:60*bci_sf,size(signal,2)-1)))); @@ -38,6 +31,6 @@ % disp(offset) %add dummy timestamp at the end to fill kinematics kin=[kin;kin(end,:)]; - kin(end,1)=size(trainingDataEEG,1)*1000*shift+windowEEG*1000; - kinPerSec=shiftingPos(kin,windowEEG,shift); + kin(end,1)=size(trainingDataEMG,1)*1000*shiftEMG+windowEMG*1000; + kinPerSec=shiftingPos(kin,windowEMG,shiftEMG); end diff --git a/usedMcode/readAllPos.m b/usedMcode/readAllPos.m index 461dd3e..ab1d427 100644 --- a/usedMcode/readAllPos.m +++ b/usedMcode/readAllPos.m @@ -1,8 +1,11 @@ -function [savePath]=readAllPos(pathToFile, subject,number,windowEMG,windowEEG,shift,maxFile,threshold,pburgOrder,minEEGFreq,maxEEGFreq,pause,noLFsamples,EMGChannels) +function [savePath]=readAllPos(pathToFile, subject,number,windowEMG,windowEEG,shiftEMG,shiftEEG,maxFile,... + threshold,pburgOrder,minEEGFreq,maxEEGFreq,pause,noLFsamples,EMGChannels) %fprintf('start read %s%i %s\n',subject,number,datestr(datetime('now'))); - savePath=strcat(pathToFile,sprintf('../matlabData/%s%i%imsWindowEMG%isWindowEEG%imsShiftFreq%ito%iPause%ipBurg%i.mat',subject,number,windowEMG*1000,windowEEG,shift*1000,minEEGFreq,maxEEGFreq,pause,pburgOrder)); - + savePath=strcat(pathToFile,... + sprintf('../matlabData/%s%i%imsWindowEMG%isWindowEEG%imsShiftEMG%imsShiftEEGFreq%ito%iPause%ipBurg%iPos.mat',... + subject,number,windowEMG*1000,windowEEG,shiftEMG*1000,shiftEEG*1000,minEEGFreq,maxEEGFreq,pause,pburgOrder)); + %only create file if it doesn't exist yet if exist(savePath, 'file') ~= 2 fprintf(strcat(savePath,' not existing; creating\n')); @@ -13,9 +16,15 @@ kin=cell([maxFile,1]); for i=1:maxFile - [sig, stat, params] = load_bcidat(strcat(pathToFile,sprintf('%s/%s_B100%i/%s_B1S00%iR0%i.dat',subject,subject,number,subject,number,i))); - tmpKin=csvread(strcat(pathToFile,sprintf('kin/%s_B1S00%iR0%i.txt',subject,number,i)),1,0); - [trainingDataEEGcell{i},trainingDataEEGlfCell{i},trainingDataEMGcell{i},kin{i}]=generateTrainingDataPos(sig,tmpKin(:,1:4),windowEMG,windowEEG,shift,params,pburgOrder,minEEGFreq,maxEEGFreq,noLFsamples,EMGChannels); + [sig, stat, params] = load_bcidat(strcat(pathToFile,... + sprintf('%s/%s_B100%i/%s_B1S00%iR0%i.dat',subject,subject,... + number,subject,number,i))); + tmpKin=csvread(strcat(pathToFile,... + sprintf('kin/%s_B1S00%iR0%i.txt',subject,number,i)),1,0); + [trainingDataEEGcell{i},trainingDataEEGlfCell{i},... + trainingDataEMGcell{i},kin{i}]=generateTrainingDataPos(sig,... + tmpKin(:,1:4),windowEMG,windowEEG,shiftEMG,shiftEEG,params,pburgOrder,... + minEEGFreq,maxEEGFreq,noLFsamples,EMGChannels); classesCell{i}=stat.StimulusCode; %fprintf('%ith file processed\n',i) end @@ -30,24 +39,39 @@ kinMat=cell2mat(kin); clear trainingDataEEGcell trainingDataEMGcell classesCell kin - classificationWithPause=classifyAccordingToEMG(size(trainingDataEEG,1), trainingDataEMG,classesMat,shift,bci_sf,threshold,pause); + classificationWithPause=classifyAccordingToEMG(trainingDataEMG,... + classesMat,shiftEMG,bci_sf,threshold,pause); clear classesMat - smoothClassification=zeros(size(classificationWithPause)); + smoothClassificationEMG=zeros(size(classificationWithPause)); for i=1:size(classificationWithPause,1) - smoothClassification(i)=round(mode(classificationWithPause(max(i-fix(bci_sf/1000),1):min(i+fix(bci_sf/1000),end)))); + smoothClassificationEMG(i)=round(mode(classificationWithPause(... + max(i-fix(bci_sf/1000),1):min(i+fix(bci_sf/1000),end)))); end clear classificationWithPause - trainingDataEEG=trainingDataEEG(smoothClassification~=-1,:,:); - trainingDataEEGlf=trainingDataEEGlf(smoothClassification~=-1,:,:); - trainingDataEMG=trainingDataEMG(smoothClassification~=-1,:); - classification=smoothClassification(smoothClassification~=-1); - kinematics=kinMat(smoothClassification~=-1,:); + smoothClassificationEEG=zeros([size(trainingDataEEG,1),1]); + for i=1:size(smoothClassificationEEG,1) + smoothClassificationEEG(i)=mode(smoothClassificationEMG(fix(... + round(shiftEEG/shiftEMG*(i-1))+1):min(fix(round(... + shiftEEG/shiftEMG*i)),max(size(smoothClassificationEMG))))); + end + + trainingDataEEG=trainingDataEEG(smoothClassificationEEG~=-1,:,:); + trainingDataEEGlf=trainingDataEEGlf(smoothClassificationEEG~=-1,:,:); + trainingDataEMG=trainingDataEMG(smoothClassificationEMG~=-1,:); + classificationEMG=smoothClassificationEMG(smoothClassificationEMG~=-1); + classificationEEG=smoothClassificationEEG(smoothClassificationEEG~=-1); + + % all created matching the size of EMGdata, if necessary downsampling + % has to be done when used + kinematics=kinMat(smoothClassificationEMG~=-1,:); + clear smoothClassification i - save(savePath,'trainingDataEEG','trainingDataEEGlf','trainingDataEMG','classification','kinematics','-v7.3'); + save(savePath,'trainingDataEEG','trainingDataEEGlf','trainingDataEMG',... + 'classificationEMG','classificationEEG','kinematics','-v7.3'); %fprintf('finished reading %s%i %s\n',subject,number,datestr(datetime('now'))); diff --git a/usedMcode/shiftingPos.m b/usedMcode/shiftingPos.m index 82dec59..dfd263b 100644 --- a/usedMcode/shiftingPos.m +++ b/usedMcode/shiftingPos.m @@ -1,10 +1,10 @@ -function [kinPerSec]=shiftingPos(kin, windowEEG, shift) - kinPerSec=zeros(fix((max(kin(:,1))-windowEEG*1000)/(shift*1000)),3); +function [kinPerSec]=shiftingPos(kin, windowEMG, shiftEMG) + kinPerSec=zeros(fix((max(kin(:,1))-windowEMG*1000)/(shiftEMG*1000)),3); for j=1:size(kinPerSec,1) - tmp=mean(kin(kin(:,1)>(j-1)*shift*1000 & kin(:,1)<=(j-1)*shift*1000+windowEEG*1000,2:4)); + tmp=mean(kin(kin(:,1)>(j-1)*shiftEMG*1000 & kin(:,1)<=(j-1)*shiftEMG*1000+windowEMG*1000,2:4)); i=j; while isnan(tmp) %interval is empty - tmp=mean(kin(kin(:,1)>i*shift*1000 & kin(:,1)<=i*shift*1000+windowEEG*1000,2:4)); + tmp=mean(kin(kin(:,1)>i*shiftEMG*1000 & kin(:,1)<=i*shiftEMG*1000+windowEMG*1000,2:4)); i=i+1; end kinPerSec(j,:)=tmp;