The case for PsychoPy

Most researchers in cognitive neuroscience have to program an experiment in a stimulus-delivery environment at some time in their career. Everyone have heard about Presentation and E-prime but much fewer have heard about PsychoPy – a free cross-platform stimulus delivery environment.

PsychoPy logo

Advantages over Presentation and E-Prime:

  1. Free and open source, meaning no license problems!
  2. Cross-platform. You can develop on Windows, Mac or Linux and run your experiment on Windows, Mac and Linux!
  3. Limitless functionality. Psychopy is a Python module, which means that it can be heavily extended with other Python modules for e.g. Matlab-level matrix operations (numpy), image and audio manipulations, advanced statistics etc.
  4. Flexible stimulus-handling: “Online” stimulus-generation, e.g. random-noise masks or moving Gaussian gratings. Also, you can define stimuli-sizes in degrees relative to the subject’s eyes, rather than having to do the trigonometry yourself.

Timing:

It compares to E-Prime and Presentation on timing. I’ve tested it on our stimulus-facilities using a photosensitive diode and comparing the timing of stimuli-appearance on the monitor to the triggers sent on a parallel port. It is reliably below 0,1 ms.

In PsychoPy you can (optionally) define durations of stimuli in monitor-frames rather than milliseconds, which ensures that timing is very well-controlled.

Disadvantages:

  1. PsychoPy’s timing of audio stimuli is still not perfect. I’ve found a ~10ms jitter in audio onset (I’ve seen much worse in some C-compiled experiments…). If you need very precise audio timing, go for something else.
  2. The serious experimenter still needs to script all of his experiment. Although a very nice graphical interface is in development, it is still not thoroughly tested and tried. In my oppinion programming in Python is much cleaner and easier than E-Prime and Presentation.

Getting started:

Good introduction to Python syntax: Learn Python in 10 minutes
Download PsychoPy at PsychoPy.org and read documentation here.

Computing half-area latency in matlab/EEGLAB

In EEG research the latency of a component (hidden in the ERP) is conventionally measured as the time when the ERP is most extreme. However, this method is not well warranted, as should be clear to anyone who have read Steven Lucks excellent chapters on ERP analysis in “An Introduction to the Event Related Potential Technique”. Instead he proposes a 50% area latency method: the latency is the point in time where the area under the curve is equal on both sides. It yields the same latencies as peak-measures on symmetrical ERPs and more correct latencies on non-symmetrical ERPs, where peak-measures fails.

A first solution

If you have your ERP as an array in MATLAB, you can compute the 50% area in this way:

To use it with EEGLAB (assuming you have your epoched EEG structure) you might implement it in this way:

How it works

Let’s have a look at how this works. Remember that we want to find the sample where the area is equal on both sides. It consists of the following steps:

  1. Transform the erp array into an array of cumulative sums, i.e. every element is the sum of all previous elements. E.g. the array erp=[1,3,4,0] becomes erp_cumsum = [1,4,8,8]. (try to imagine this ERP to ease your understanding). Each element then corresponds to the area that precede it in the ERP. The last element corresponds the total area of the ERP.
  2. Now the job is simply to find the midmost value in erp_cumsum, i.e. the value that separates the area into two equal parts. To do this, we subtract the half area from every element in erp_cumsum. Continuing the example from above, sum(erp)/2 = 4. The subtraction giveserp_reduced = erp_cumsum-sum(erp)/2 = [-3,0,5,5].
  3. The midmost value in erp_cumsum is now the element closest to zero (we can already see that this is sample number 2). Finding it in matlab is done by taking the absolute value, which gives erp_abs = abs(erp_reduced) = [3,0,5,5], and finding the index of the minimum value. This is simply [~,index] = min(erp_abs).

Correction: defining the direction of the component

This method is not fool-proof but things can be done to prevent some errors. Consider a positive ERP (e.g. P300) that starts and ends at negative values, .e.g erp = [-3,0,5,2,-1,-4]. This is a very nice peak and visually you would consider it all to be part of the P300. However, using cumulative sums, the start and end would be “negative area”, which doesn’t make sense when we’re trying to find a latency.

The solution starts off by defining whether you’re looking for the latency of a positive or a negative peak. Then you adjust your ERP to only contain positive or negative values by raising or lowering the whole ERP so that all elements have the same sign. In this example, we’re looking for a positive component and then we raise the ERP by 4. If we were looking for a negative component, we’d lower the ERP by 5. So:

Extending it to multiple ERPs

Say you’re dealing with multiple ERPs and want to calculate 50% latencies on all of them. You want to put in e.g. 50 ERPs and get 50 latencies as output. This is practical if you want to calculate latencies on individual trials to investigate variability or on multiple subjects at the same time. This is where Matlab’s excellent matrix operations get to work. I don’t want to go into details on this one. I’m just going to give you the code for positive-going components in single-trial ERPs: