Olfaction#

Summary: In this tutorial, we will implement a simple controller for odor-guided taxis.

In the previous tutorial, we covered how one can simulate the visual experience of the fly simulation. In addition to vision, we also made it possible for our model to detect odors emitted by objects in the simulation environment. The olfactory system in Drosophila consists of specialized olfactory sensory neurons (OSNs) located in the antennae and maxillary palps. These detect specific odorant molecules and convey this information to the brain’s antennal lobe, where their signals are further processed. This is shown in the figure below (left, source: Martin et al, 2013) We emulated peripheral olfaction by attaching virtual odor sensors to the antennae and maxillary palps of our biomechanical model, as shown in the figure (right). The user has the option of configuring additional sensors at more precise locations on these olfactory organs. These virtual sensors can detect odor intensities across a multi-dimensional space that can be thought of as representing, for example, the concentrations of monomolecular chemicals sensed by OSNs in the antennae, or the intensities of composite odors co-activating numerous projection neurons in the antennal lobe.

https://github.com/NeLy-EPFL/_media/blob/main/flygym/olfaction.png?raw=true

Odor arena#

To demonstrate odor sensing, let’s create an environment with one attractive odor source and two aversive odor sources. The dimension of this odor space is 2 (attractive, aversive) despite the number of odor sources being 3. The odor sources share a peak intensity of 1. We will color the attractive odor source orange and the aversive odor sources blue.

 1 import numpy as np
 2
 3 # Odor source: array of shape (num_odor_sources, 3) - xyz coords of odor sources
 4 odor_source = np.array([[24, 0, 1.5], [8, -4, 1.5], [16, 4, 1.5]])
 5
 6 # Peak intensities: array of shape (num_odor_sources, odor_dimesions)
 7 # For each odor source, if the intensity is (x, 0) then the odor is in the 1st dimension
 8 # (in this case attractive). If it's (0, x) then it's in the 2nd dimension (in this case
 9 # aversive)
10 peak_intensity = np.array([[1, 0], [0, 1], [0, 1]])
11
12 # Marker colors: array of shape (num_odor_sources, 4) - RGBA values for each marker,
13 # normalized to [0, 1]
14 marker_colors = [[255, 127, 14], [31, 119, 180], [31, 119, 180]]
15 marker_colors = np.array([[*np.array(color) / 255, 1] for color in marker_colors])
16
17 odor_dimesions = len(peak_intensity[0])

Let’s create the arena using these parameters. The detailed documentation of the OdorArena class can be found in the API reference. Its implementation is beyond the scope of this tutorial but can be found here.

1 from flygym.mujoco.arena import OdorArena
2
3 arena = OdorArena(
4     odor_source=odor_source,
5     peak_intensity=peak_intensity,
6     diffuse_func=lambda x: x**-2,
7     marker_colors=marker_colors,
8     marker_size=0.3,
9 )

Let’s place our fly in the arena. As before, we will run a few iterations to allow it to stand on the ground in a stable manner.

 1 import matplotlib.pyplot as plt
 2 from flygym.mujoco import Parameters
 3 from flygym.mujoco.examples.turning_controller import HybridTurningNMF
 4
 5
 6 contact_sensor_placements = [
 7     f"{leg}{segment}"
 8     for leg in ["LF", "LM", "LH", "RF", "RM", "RH"]
 9     for segment in ["Tibia", "Tarsus1", "Tarsus2", "Tarsus3", "Tarsus4", "Tarsus5"]
10 ]
11 sim_params = Parameters(
12     timestep=1e-4,
13     render_mode="saved",
14     render_playspeed=0.5,
15     render_window_size=(800, 608),
16     enable_olfaction=True,
17     enable_adhesion=True,
18     draw_adhesion=False,
19     render_camera="birdeye_cam",
20 )
21 sim = HybridTurningNMF(
22     sim_params=sim_params,
23     arena=arena,
24     spawn_pos=(0, 0, 0.2),
25     contact_sensor_placements=contact_sensor_placements,
26 )
27 for i in range(500):
28     sim.step(np.zeros(2))
29     sim.render()
30 fig, ax = plt.subplots(1, 1, figsize=(5, 4), tight_layout=True)
31 ax.imshow(sim._frames[-1])
32 ax.axis("off")
33 fig.savefig("./outputs/olfaction_env.png")
https://github.com/NeLy-EPFL/_media/blob/main/flygym/olfaction_env.png?raw=true

Controller for odor taxis#

Let’s design a simple hand-tuned controller for odor-guided taxis. We start by calculating the left-right asymmetry of the odor intensity \(I\) for each odor \(o\):

\[\Delta I_o = \frac{I_\text{left,o} - I_\text{right,o}}{(I_\text{left,o} + I_\text{right,o}) / 2}\]

Then, we multiply \(\Delta I_o\) by a gain \(\gamma_o\) for each odor dimension and take the sum \(s\). Attractive and aversive odors will have different signs in their gains.

\[s = \sum_{o} \gamma_o \Delta I_o\]

We transform \(s\) nonlinearly to avoid overly drastic turns when the asymmetry is subtle and to crop it within the range [0, 1). This gives us a turning bias \(b\):

\[b = \tanh(s^2)\]

Finally, we modulate the descending signal \(\delta\) based on \(b\) and the sign of \(s\):

\[\begin{split}\delta_\text{left} = \begin{cases} \delta_\text{max} & \text{if } s>0\\ \delta_\text{max} - b(\delta_\text{max} - \delta_\text{min}) & \text{otherwise} \end{cases} \qquad \delta_\text{right} = \begin{cases} \delta_\text{max} - b(\delta_\text{max} - \delta_\text{min}) & \text{if } s>0\\ \delta_\text{max} & \text{otherwise} \end{cases}\end{split}\]

where, \(\delta_\text{min}\), \(\delta_\text{max}\) define the range of the descending signal. Here, we will use the following parameters:

  • \(\gamma_\text{attractive} = -500\) (negative ipsilateral gain leads to positive taxis)

  • \(\gamma_\text{aversive} = 80\) (positive ipsilateral gain leads to negative taxis)

  • \(\delta_\text{min} = 0.2\)

  • \(\delta_\text{max} = 1\)

As before, we will recalculate the steering signal every 0.05 seconds. Let’s implement this in Python:

 1 from tqdm import trange
 2
 3 attractive_gain = -500
 4 aversive_gain = 80
 5 decision_interval = 0.05
 6 run_time = 5
 7 num_decision_steps = int(run_time / decision_interval)
 8 physics_steps_per_decision_step = int(decision_interval / sim_params.timestep)
 9
10 obs_hist = []
11 odor_history = []
12 obs, _ = sim.reset()
13 for i in trange(num_decision_steps):
14     attractive_intensities = np.average(
15         obs["odor_intensity"][0, :].reshape(2, 2), axis=0, weights=[9, 1]
16     )
17     aversive_intensities = np.average(
18         obs["odor_intensity"][1, :].reshape(2, 2), axis=0, weights=[10, 0]
19     )
20     attractive_bias = (
21         attractive_gain
22         * (attractive_intensities[0] - attractive_intensities[1])
23         / attractive_intensities.mean()
24     )
25     aversive_bias = (
26         aversive_gain
27         * (aversive_intensities[0] - aversive_intensities[1])
28         / aversive_intensities.mean()
29     )
30     effective_bias = aversive_bias + attractive_bias
31     effective_bias_norm = np.tanh(effective_bias**2) * np.sign(effective_bias)
32     assert np.sign(effective_bias_norm) == np.sign(effective_bias)
33
34     control_signal = np.ones((2,))
35     side_to_modulate = int(effective_bias_norm > 0)
36     modulation_amount = np.abs(effective_bias_norm) * 0.8
37     control_signal[side_to_modulate] -= modulation_amount
38
39     for j in range(physics_steps_per_decision_step):
40         obs, _, _, _, _ = sim.step(control_signal)
41         rendered_img = sim.render()
42         if rendered_img is not None:
43             # record odor intensity too for video
44             odor_history.append(obs["odor_intensity"])
45         obs_hist.append(obs)
46
47     # Stop when the fly is within 2mm of the attractive odor source
48     if np.linalg.norm(obs["fly"][0, :2] - odor_source[0, :2]) < 2:
49         break
77%|███████▋  | 77/100 [01:48<00:32,  1.41s/it]

We can visualize the fly trajectory:

 1 fly_pos_hist = np.array([obs["fly"][0, :2] for obs in obs_hist])
 2 fig, ax = plt.subplots(1, 1, figsize=(5, 4), tight_layout=True)
 3 ax.scatter(
 4     [odor_source[0, 0]],
 5     [odor_source[0, 1]],
 6     marker="o",
 7     color="tab:orange",
 8     s=50,
 9     label="Attractive",
10 )
11 ax.scatter(
12     [odor_source[1, 0]],
13     [odor_source[1, 1]],
14     marker="o",
15     color="tab:blue",
16     s=50,
17     label="Aversive",
18 )
19 ax.scatter([odor_source[2, 0]], [odor_source[2, 1]], marker="o", color="tab:blue", s=50)
20 ax.plot(fly_pos_hist[:, 0], fly_pos_hist[:, 1], color="k", label="Fly trajectory")
21 ax.set_aspect("equal")
22 ax.set_xlim(-1, 25)
23 ax.set_ylim(-5, 5)
24 ax.set_xlabel("x (mm)")
25 ax.set_ylabel("y (mm)")
26 ax.legend(ncols=3, loc="lower center", bbox_to_anchor=(0.5, -0.6))
27 fig.savefig("./outputs/odor_taxis_trajectory.png")
https://github.com/NeLy-EPFL/_media/blob/main/flygym/odor_taxis_trajectory.png?raw=true

We can also generate the video:

1 sim.save_video("./outputs/odor_taxis.mp4")