Automazione funzioni

di il
2 risposte

Automazione funzioni

Ciao a tutti,
premetto che sono nuova su Matlab ma devo usarlo per lavoro e non ho molta idea di cosa sto scrivendo come codice.
In pratica devo scrivere un codice che mi permetta di dividere il segnale che sto registrando con degli elettrodi in base al numero di elettrodi stesso. Gli elettrodi servono a registrare attività cerebrale e l'obiettivo è riuscire a isolare gli spikes mandati dai neuroni che rispondono ad un dato stimolo acustico. Vi riporto sotto il codice che ho scritto estraendo script già usati dal mio capo. Tutto funziona e riesco sia a calcolare gli spikes per ogni elettrodo (definito come canale) che a generare i grafici che ci servono a capire come posizionare gli elettrodi. L'unico problema è che per eseguire i calcoli su ogni canale io ho riportato il nome del file di registrazione per ogni canale e invece vorrei riuscire ad automatizzare il tutto in maniera tale da non dover eseguire n volte le stesse operazioni su n canali ma eseguirne solo un blocco che mi funzioni su tutti gli n canali.
Qualsiasi suggerimento è ben accetto!
Sotto ho riportato solo i codici per due canali ma il concetto credo si sia capito...ho fatto un copia incolla per tutti i canali.
La piattaforma che registra l'attività attraverso i 16 canali genera 16 file che si chiamano 100_CH0.continuous....fino a 100_CH15.continuous quindi io ho rinominato ogni volta le variabili che genero (data, info, suprasigma etc) con l'indice del numero del canale..

come faccio ad automatizzare tutto?

GRAZIE MILLE!

do_plot = 1;
%channel 0
p = inputParser;
Validchar = @(x) assert(ischar(x),'oefile');
validScalarPosNum = @(x) assert(isnumeric(x) && isscalar(x) && (x > 0),'1');
LogicalValue = @(x) assert(islogical(x),'false');
addRequired(p,'oefile',Validchar);
addOptional(p,'sigmatresh',5,validScalarPosNum);
addOptional(p,'bandpass_lower',150,validScalarPosNum);
addOptional(p,'bandpass_upper',3000,validScalarPosNum);
addParameter(p,'GraphicalOutput',false,LogicalValue);
parse(p,'100_CH0.continuous');

[data0,~,info0] = load_open_ephys_data('100_CH0.continuous');

sigmatresh = p.Results.sigmatresh; bandpass_lower = p.Results.bandpass_lower; bandpass_upper = p.Results.bandpass_upper;
bp_filtered = bandpass(data0, [bandpass_lower bandpass_upper],info0.header.sampleRate);
suprasigma0 = find(bp_filtered<-median(abs(bp_filtered)/0.6745)*sigmatresh);
suprasigma0(diff(suprasigma0)==1) = NaN;
suprasigma0 = suprasigma0(~isnan(suprasigma0));

for i = 1:length(suprasigma0)
sigmatrixindex0(i,:) = suprasigma0(i)-20:suprasigma0(i)+20;
end

sigmatrix0 = bp_filtered(sigmatrixindex0);
[~,linindex] = min(sigmatrix0,[],2,'omitnan');

for i = 1:length(linindex)
sigmatrixindex_aligned0(i,:) = sigmatrixindex0(linindex(i))-225:sigmatrixindex0(linindex(i))+225;
end
sigmatrix_aligned0 = bp_filtered(sigmatrixindex_aligned0);

if do_plot==1
figure
set(gcf,'color','w');
subplot(221)
histogram(suprasigma0,'LineWidth',1.2)
title('Histogram spikes events CH0')
ylabel('number of spikes')
xlabel('time (ms)')
xticklabels(num2cell((xticks/info0.header.sampleRate)*1000))
subplot(222)
plot(mean(sigmatrix_aligned0),'k','LineWidth',1.2)
title('Meaned spike trace CH0')
ylabel('mV')
xlabel('time (ms)')
xticklabels(num2cell((xticks/info0.header.sampleRate)*1000))

subplot(2,2,[3 4])
plot(bp_filtered,'-*','Color','k','LineWidth',0.1,'MarkerEdgeColor','red','MarkerIndices',suprasigma0)
title('Bandpass filtered signal with highlighted spikes CH0')
ylabel('mV')
xlabel('time (s)')
xticklabels(xticks/info0.header.sampleRate)
end

spikecount0 = length(suprasigma0);
spikerate0 = spikecount0/(length(bp_filtered)/info0.header.sampleRate);


%channel 1

p = inputParser;
Validchar = @(x) assert(ischar(x),'oefile');
validScalarPosNum = @(x) assert(isnumeric(x) && isscalar(x) && (x > 0),'1');
LogicalValue = @(x) assert(islogical(x),'false');
addRequired(p,'oefile',Validchar);
addOptional(p,'sigmatresh',5,validScalarPosNum);
addOptional(p,'bandpass_lower',150,validScalarPosNum);
addOptional(p,'bandpass_upper',3000,validScalarPosNum);
addParameter(p,'GraphicalOutput',false,LogicalValue);
parse(p,'100_CH1.continuous');

[data1,~,info1] = load_open_ephys_data('100_CH1.continuous');

sigmatresh = p.Results.sigmatresh; bandpass_lower = p.Results.bandpass_lower; bandpass_upper = p.Results.bandpass_upper;
bp_filtered = bandpass(data1, [bandpass_lower bandpass_upper],info1.header.sampleRate);
suprasigma1 = find(bp_filtered<-median(abs(bp_filtered)/0.6745)*sigmatresh);
suprasigma1(diff(suprasigma1)==1) = NaN;
suprasigma1 = suprasigma1(~isnan(suprasigma1));

for i = 1:length(suprasigma1)
sigmatrixindex1(i,:) = suprasigma1(i)-20:suprasigma1(i)+20;
end

sigmatrix1 = bp_filtered(sigmatrixindex1);
[~,linindex] = min(sigmatrix1,[],2,'omitnan');

for i = 1:length(linindex)
sigmatrixindex_aligned1(i,:) = sigmatrixindex1(linindex(i))-225:sigmatrixindex1(linindex(i))+225;
end
sigmatrix_aligned1 = bp_filtered(sigmatrixindex_aligned1);

if do_plot==1
figure
subplot(221)
set(gcf,'color','w');
histogram(suprasigma1,'LineWidth',1.2)
title('Histogram spikes events CH1')
ylabel('number of spikes')
xlabel('time (ms)')
xticklabels(num2cell((xticks/info1.header.sampleRate)*1000))
subplot(222)
plot(mean(sigmatrix_aligned1),'k','LineWidth',1.2)
title('Meaned spike trace CH1')
ylabel('mV')
xlabel('time (ms)')
xticklabels(num2cell((xticks/info1.header.sampleRate)*1000))

subplot(2,2,[3 4])
plot(bp_filtered,'-*','Color','k','LineWidth',0.1,'MarkerEdgeColor','red','MarkerIndices',suprasigma1)
title('Bandpass filtered signal with highlighted spikes CH1')
ylabel('mV')
xlabel('time (s)')
xticklabels(xticks/info1.header.sampleRate)
end
spikecount1 = length(suprasigma1);
spikerate1 = spikecount1/(length(bp_filtered)/info1.header.sampleRate);

2 Risposte

  • Re: Automazione funzioni

    Considerando che la struttura del nome dei files è sempre la stessa, la cosa più semplice che puoi fare è:

    [*] racchiudere il codice che corrisponde ad un caso all'interno di un ciclo "for"
    [*] generare il nome del file di input come stringa in modo automatico
    [*] sostituire, nel blocco di codice, tutte le occorrenze del nome del file (che al momento scrivi direttamente) con la stringa che hai generato

    La struttura del nome dei files si compone di tre parti: la prima e la terza sono costanti, la seconda è un numero intero progressivo (XXX)

    f_name = 100_CH
    XXX
    .continuous

    Per costruire la stringa con il nome dei file puoi definire due stringhe che contengono la prima e la terza parte del nome del file; la seconda parte (il numero progressivo) puoi crearla con la funzione "num2str" che converte un numero in una stringa
    Per creare la stringa completa basta concatenare le tre parti come nell'esempio che segue che crea semplicemente i nomi dei 16 file.
    
    % Prima parte del nome del file di input
    base_name_1='100_CH'
    % Terza parte del nome del file di input
    base_name_3='.continuous'
    
    for i=0:15
       f_name=[base_name_1 num2str(i) base_name_3]
       
       %
       % INSERISCI QUI IL BLOCCO DI CODICE
       %
       
    end
    
    A questo punto sarà sufficiente sostituire le occorrenze del nome del file nel blocco di codice con la stringa.
    Di seguito un estratto del tuo codice, modificato:
    
    %parse(p,'100_CH1.continuous');
    parse(p,f_name);
    
    %[data1,~,info1] = load_open_ephys_data('100_CH1.continuous');
    [data1,~,info1] = load_open_ephys_data(f_name);
    
    
    Se racchiudi il blocco di codice all'interno di un loop, i 16 casi verranno eseguiti in sequenza ed il codice genererà 16 figure (ricorda di inserire la chiamata alla funzione "figure" nel blocco).

    Qualche considerazione sul codice che hai pubblicato:
    [*] ammesso che per ogni blocco tu abbia inserito tutto il codice, non c'è un punto nel quale salvi i risultati ottenuti (le figure o le variabili) da usare per successive analisi (ammesso che ne vengano fatte)
    [*] al momento usi variabili diverse in ogni blocco (es: "suprasigma0" e "suprasigma1", "sigmatrix0" e "sigmatrix1")
    se utilizzerai un solo blocco di codice all'intero di un loop, non potrai avere, nel workspace i valori delle variabili di tutti i casi (ad esempio le 16 sigmatrix0, sigmatrix1, ... o qualunque altra variabile)
    Nel caso ti servissero, devi salvare le variabili in vettori o matrici per poterle usare al termine del loop.

    Vista la ripetitività del codice, potresti anche creare un'interfaccia grafica (con il tool "guide" o con "app designer") che potrebbe contenere, ad esempio:
    [*] un listbox per la selezione dei casi / file di input (questo ti permetterebbe di eseguire un caso specifico)
    [*] un axes nel quale plottare il caso selezionato
    [*] pushbuttons per avviare l'analisi di un caso o di tutti (in questo caso le figure non verrebbero plottate nell'interfaccia ma in "figure" separate)
    [*] dei pushbutton per selezionare file di input da cartelle diverse (ammesso che ti serva)
    [*] dei pushbuttons per salvare i risultati dell'analisi (figure, dati)
    [*] dei pushbuttons per aprite casi elaborati e salvati in precedenza
    e tante altre cose.
  • Re: Automazione funzioni

    Ciao,
    ho provato il tuo codice e funziona perfettamente!
    Grazie infinite davvero!
    Gentilissimo
Devi accedere o registrarti per scrivere nel forum
2 risposte