Skeleton code for automatic fits


Introduction

One way to save yourself a lot of time is to write analysis code that can automatically perform a first-pass analysis of an experiment, by 1) generating fits and calculating index values that are commonly used, and 2) generating simple plots so that you can examine the data. When I was a postdoc and had things working smoothly, I would run these automatic functions while I was cleaning up the experiment, and by the next morning the data would be ready to interpret.

When to run this type of analysis

I run these analyses after all the primary data collection is finished. That is
  • the log is complete
  • any little errors that occurred during data collection have been addressed
  • any directory that has improperly or poorly collected data (that is, a cell died or the image went blank) has been "neutered" by changing the reference.txt filename to reference0.txt or something similar.
  • spikes have been extracted, or fluorescence has been calculated from images

A code example of a function that processes an experiment

This is analyzetftraining.m, a file that analyzes the impact of motion training at one speed on direction selectivity at many speeds:

function [cells,cellnames]=analyzetftraining(thedir, redoresponses, redoanalysis, plotit, saveit)
% ANALYZESCRAMBLED - Perform fitting analysis for a tftraining training experiment
%
%  [CELLS,CELLNAMES]=ANALYZETFTRAINING(THEDIR, REDORESPONSES, REDOANALYSIS,...
%    PLOTIT, SAVEIT)
%
%  Analyses data from the experiment in directory THEDIR.
%  Inputs:
%  DORESPONSES - If is 1, then the raw responses are re-extracted
%  REDOANALSYS - If is 1, then fits are performed
%  PLOTIT - If is 1, then a plot is made for each extraction or fit
%  SAVEIT - If 1, the extractions/fit results are saved to the experiment.mat file

  % STEP 1: READ IN THE CELLS

ds = dirstruct(thedir);

[cells,cellnames] = load2celllist(getexperimentfile(ds),'cell*','-mat');

allcellnames = cellnames;

  % STEP 2: ADD INFORMATION THAT IDENTIFIES THE STIMULI IN EACH DIRECTORY

[dummy,cells] = add_testdir_info(ds,cells); % add the directory information

  % STEP 3: If we are importing spikes from Plexon, we need our own records of unit quality

%cellinfo = read_unitquality(ds); % not needed since we're not importing from plexon
%[cells,cellnames]=filter_by_quality(ds,cells,cellnames,cellinfo); % limit by cell quality

  % STEP 4: Limit the cells we will analyze to those that are numbered between 50 and 75
  %   for Neil's experiment, we want only the spikes sorted with the VH lab tool

[cells,cellnames] = filter_by_index(cells,cellnames,50,75); % only get VH lab sorted spikes

  % STEP 5: Give the user a chance to see the names of the cells that are going to be analyzed
  %          this gives the user a chance to abort before hours of analysis run

disp(['Cells we have now:']);
cellnames,

pause(5);

% STEP 6: Identify the type of training that was used

training_assoc = read_trainingtype(ds,'ErrorIfNoTrainingType',1,'ErrorIfNoTrainingAngle',1,...
                'ErrorIfNoTF',1);

 % STEP 7: extract responses -- use spike times and stimulus times to determine stimulus responses

    % get the list of all the stimulus types we know
[dirtestnamelist,velocitytestnamelist,tftrainingtestnamelist] = vhdirectiontrainingtypes(1);

if redoresponses, 
    for i=1:length(cells),
        disp(['Now analyzing cell ... ' cellnames{i} ', (' int2str(i) ' of ' int2str(length(cells)) ')']);
        for j=1:length(training_assoc),
                 % add training info to each cell's information in the database
            cells{i} = associate(cells{i},training_assoc(j));
        end;
        % identify the time of each recording
        cells{i}=extractstimdirectorytimes(ds,cells{i},...
               'ErrorIfEmpty',1,'EarlyMorningCutOffTime',7);
        % Now do the extraction for all stimulus types we know: See extraction notes below
        cells{i}=performsingleunitgrating(ds,cells{i},cellnames{i},...
            dirtestnamelist,'angle',plotit);
        cells{i}=performsingleunitgrating(ds,cells{i},cellnames{i},...
              velocitytestnamelist,'tFrequency',plotit);
        cells{i}=performsingleunitgrating(ds,cells{i},cellnames{i},...
              tftrainingtestnamelist,'stimnumber',plotit);
    end;

    if saveit,
        saveexpvar(ds,cells,cellnames);
    else,
        disp(['Not saving per user request.']);
    end;

end;

if redoanalysis,
    for i=1:length(cells),
        disp(['Now fitting cell ... ' cellnames{i} ...
            ', (' int2str(i) ' of ' int2str(length(cells)) ')']);
        % perform orientation and direction fitting; see fitting notes below
        cells{i} = performspperiodicgenericanalysis(ds,cells{i},cellnames{i},...
                plotit,dirtestnamelist,'p.angle',...
                'otanalysis_compute','',3);
        
    end;
    if saveit,
        saveexpvar(ds,cells,cellnames);
    else,
        disp(['Not saving per user request.']);
    end;
end;


Extraction notes:

[Notes here on performsingleunitgrating]

Fitting notes:

[Notes here on performspperiodicgenericanalysis]


Comments