Electrical Stimuli

The Stimulus object defines a common interface for all electrical stimuli. It provides a 2-D data array with labeled axes, where rows denote electrodes and columns denote points in time. Stimuli can be assigned to electrodes of a ProsthesisSystem object, who will deliver them to the retina.

A stimulus can be created from a variety of source types, such as the following:

  • Scalar value: interpreted as the current amplitude delivered to a single electrode (no time component).
  • NumPy array:
    • Nx1 array: interpreted as N current amplitudes delivered to N electrodes (no time component).
    • NxM array: interpreted as N electrodes each receiving M current amplitudes in time.
  • TimeSeries: interpreted as the stimulus in time for a single electrode (e.g., PulseTrain).

In addition, you can also pass a collection of source types (e.g., list, dictionary).

See the Stimulus API documentation for a full list.

Note

Depending on the source type, a stimulus might have a time component or not.

Single-electrode stimuli

The easiest way to create a stimulus is to specify the current amplitude (uA) to be delivered to an electrode:

In [1]: from pulse2percept.stimuli import Stimulus

# Stimulate an unnamed electrode with -14uA:
In [2]: Stimulus(-14)
Out[2]: 
Stimulus(data=[[-14.]], electrodes=[0], metadata=None, 
         shape=(1, 1), time=None)

You can also specify the name of the electrode to be stimulated:

# Stimulate Electrode 'B7' with -14uA:
In [3]: Stimulus(-14, electrodes='B7')
Out[3]: 
Stimulus(data=[[-14.]], electrodes=['B7'], metadata=None, 
         shape=(1, 1), time=None)

By default, this stimulus will not have a time component (stim.time is None). Some models, such as ScoreboardModel, cannot handle stimuli in time.

To create stimuli in time, you can pass a TimeSeries object, such as a BiphasicPulse or a PulseTrain:

# Stimulate Electrode 'A001' with a cathodic-first 20Hz pulse train
# with 10uA amplitude, lasting for 0.5s, sampled at 0.1ms:
In [4]: from pulse2percept.stimuli import PulseTrain

In [5]: pt = PulseTrain(0.0001, freq=20, amp=10, dur=0.5)

In [6]: stim = Stimulus(pt)

In [7]: stim
Out[7]: 
Stimulus(data=[[-10. -10. ...   0.   0.]], electrodes=[0], 
         metadata=None, shape=(1, 5000), 
         time=[0.e+00 1.e-04 ... 5.e-01 5.e-01])

# This stimulus has a time component:
In [8]: stim.time
Out[8]: array([0.e+00, 1.e-04, ..., 5.e-01, 5.e-01])

You can specify not only the name of the electrode but also the time steps to be used:

# Stimulate Electrode 'C7' with int time steps:
In [9]: Stimulus(pt, electrodes='C7', time=np.arange(pt.shape[-1]))
Out[9]: 
Stimulus(data=[[-10. -10. ...   0.   0.]], 
         electrodes=['C7'], metadata=None, shape=(1, 5000), 
         time=[   0    1 ... 4998 4999])

Creating multi-electrode stimuli

Stimuli can also be created from a list or dictionary of source types:

# Stimulate three unnamed electrodes with -2uA, 14uA, and -100uA,
# respectively:
In [10]: Stimulus([-2, 14, -100])
Out[10]: 
Stimulus(data=[[  -2.], [  14.], [-100.]], 
         electrodes=[0 1 2], metadata=None, shape=(3, 1), 
         time=None)

Electrode names can be passed in a list:

In [11]: Stimulus([-2, 14, -100], electrodes=['A1', 'B1', 'C1'])
Out[11]: 
Stimulus(data=[[  -2.], [  14.], [-100.]], 
         electrodes=['A1' 'B1' 'C1'], metadata=None, 
         shape=(3, 1), time=None)

Alternatively, stimuli can be created from a dictionary:

# Equivalent to the previous one:
In [12]: Stimulus({'A1': -2, 'B1': 14, 'C1': -100})
Out[12]: 
Stimulus(data=[[  -2.], [  14.], [-100.]], 
         electrodes=['A1' 'B1' 'C1'], metadata=None, 
         shape=(3, 1), time=None)

The same is true for a dictionary of pulse trains:

# Sending the same pulse train to three specific electrodes:
In [13]: Stimulus({'A1': pt, 'B1': pt, 'C1': pt})
Out[13]: 
Stimulus(data=<(3, 5000) np.ndarray>, 
         electrodes=['A1' 'B1' 'C1'], metadata=None, 
         shape=(3, 5000), 
         time=[0.e+00 1.e-04 ... 5.e-01 5.e-01])

Assigning new coordinates to an existing stimulus

You can change the coordinates of an existing Stimulus object, but retain all its data, as follows:

# Say you have a Stimulus object with unlabeled axes:
In [14]: stim = Stimulus(np.ones((2, 5)))

In [15]: stim
Out[15]: 
Stimulus(data=[[1. 1. ... 1. 1.], [1. 1. ... 1. 1.]], 
         electrodes=[0 1], metadata=None, shape=(2, 5), 
         time=[0 1 2 3 4])

# You can create a new object from it with named electrodes:
In [16]: Stimulus(stim, electrodes=['A1', 'F10'])
Out[16]: 
Stimulus(data=[[1. 1. ... 1. 1.], [1. 1. ... 1. 1.]], 
         electrodes=['A1' 'F10'], metadata=None, 
         shape=(2, 5), time=[0 1 2 3 4])

# Same goes for time points:
In [17]: Stimulus(stim, time=[0, 0.1, 0.2, 0.3, 0.4])
Out[17]: 
Stimulus(data=[[1. 1. ... 1. 1.], [1. 1. ... 1. 1.]], 
         electrodes=[0 1], metadata=None, shape=(2, 5), 
         time=[0.  0.1 0.2 0.3 0.4])

Compressing a stimulus

The compress method automatically compresses the data in two ways:

  • Removes electrodes with all-zero activation.
  • Retains only the time points at which the stimulus changes.

For example, only the signal edges of a pulse train are saved. That is, rather than saving the current amplitude at every 0.1ms time step, only the non-redundant values are retained. This drastically reduces the memory footprint of the stimulus. You can convince yourself of that by inspecting the size of a Stimulus object before and after compression:

# An uncompressed stimulus:
In [18]: stim = Stimulus(PulseTrain(0.0001, freq=10), compress=False)

In [19]: stim
Out[19]: 
Stimulus(data=[[-20. -20. ...   0.   0.]], electrodes=[0], 
         metadata=None, shape=(1, 5000), 
         time=[0.e+00 1.e-04 ... 5.e-01 5.e-01])

# Now compress the data:
In [20]: stim.compress()

# Notice how the stimulus shape and time axis have changed:
In [21]: stim
Out[21]: 
Stimulus(data=[[-20. -20. ...   0.   0.]], electrodes=[0], 
         metadata=None, shape=(1, 40), 
         time=[0.00e+00 3.00e-04 ... 4.01e-01 5.00e-01])

Interpolating stimulus values

The Stimulus object can also interpolate stimulus values at time points that are not explicitly provided. Interpolation is done automatically by providing an index or slice into the 2-D data array. The first dimension addresses the electrode (and cannot be interpolated), and the second dimension addresses time (which can be interpolated):

# A single-electrode ramp stimulus:
In [22]: stim = Stimulus(np.arange(10).reshape((1, -1)))

In [23]: stim
Out[23]: 
Stimulus(data=[[0. 1. ... 8. 9.]], electrodes=[0], 
         metadata=None, shape=(1, 10), time=[0 1 ... 8 9])

# Retrieve stimulus at t=3:
In [24]: stim[0, 3]
Out[24]: 3.0

# Time point 3.45 is not in the data provided above, but can be
# interpolated as follows:
In [25]: stim[0, 3.45]
Out[25]: 3.45

# This also works for multiple time points:
In [26]: stim[0, [3.45, 6.78]]
Out[26]: 3.45

# Extrapolating is disabled by default, but you can enable it:
In [27]: stim = Stimulus(np.arange(10).reshape((1, -1)), extrapolate=True)

In [28]: stim[0, 123.45]
Out[28]: 123.45

For a multi-electrode stimulus, you can access stimulus values at time t for any or all electrodes:

# Multi-electrode stimulus
In [29]: stim = Stimulus(np.arange(100).reshape((5, 20)))

In [30]: stim
Out[30]: 
Stimulus(data=<(5, 20) np.ndarray>, electrodes=[0 1 2 3 4], 
         metadata=None, shape=(5, 20), 
         time=[ 0  1 ... 18 19])

# Interpolate t=4.5 for all electrodes:
In [31]: stim[:, 4.5]
Out[31]: 
array([[ 4.5],
       [24.5],
       [44.5],
       [64.5],
       [84.5]])

You can choose different interpolation methods, as long as scipy.interpolate.interp1d accepts them. For example, the ‘nearest’ method will return the value of the nearest data point:

# A single-electrode ramp stimulus:
In [32]: stim = Stimulus(np.arange(10).reshape((1, -1)), interp_method='nearest',
   ....:                 extrapolate=True)
   ....: 

# Interpolate:
In [33]: stim[0, 3.45]
Out[33]: 3.0

# Outside the data range:
In [34]: stim[0, 12.2]
Out[34]: 9.0